Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Rate neuron framework #745

Merged
merged 47 commits into from
Aug 9, 2017
Merged

Rate neuron framework #745

merged 47 commits into from
Aug 9, 2017

Conversation

janhahne
Copy link
Contributor

@janhahne janhahne commented Jun 6, 2017

This PR adds rate models to NEST as described in

Hahne, J., Dahmen, D., Schuecker, J., Frommer, A.,
Bolten, M., Helias, M. and Diesmann, M. (2017).
Integration of Continuous-Time Dynamics in a
Spiking Neural Network Simulator.
Front. Neuroinform. 11:34. doi: 10.3389/fninf.2017.00034

Overview:

  • added rate model base classes for rate neurons with input noise (rate_neuron_ipn) and outout noise (rate_neuron_opn)
  • added rate connections (rate_connection for connections without delay, delay_rate_connection for connections with delay and diffussion_connection for connection between siegert neurons)
  • added several template based rate neuron models (lin_rate, tanh_rate, thresholdlin_rate)
  • added siegert_neuron model
  • added SecondaryEvents for different rate connection
  • adjusted some existing unittests to work with or exclude rate neurons
  • added sli and python unittests (functionality, mpi, connections) for the new framework
  • added two python examples

We suggest @heplesser as a reviewer

@heplesser heplesser requested review from heplesser and mschmidt87 June 6, 2017 12:57
@heplesser heplesser added ZC: Kernel DO NOT USE THIS LABEL I: Internal API Changes were introduced in basic internal workings of the simulator that developers need to know ZP: PR Created DO NOT USE THIS LABEL S: High Should be handled next T: Enhancement New functionality, model or documentation labels Jun 6, 2017
@janhahne janhahne force-pushed the rate-neuron-framework branch from 709a268 to 728436c Compare June 6, 2017 13:50
Copy link

@mschmidt87 mschmidt87 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Very impressive work.
So far, I have added some comment and requests, which are mostly referring to documentation issues.

between neurons of type siegert_neuron. The connection type is identical to
type rate_connection for instantaneous rate connections except for the two
parameters "drift_factor" and "diffusion_factor" substituting the
parameter "weight".

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Perhaps you could add one more sentence explaining the meaning of drift and diffusion in the context of the siegert neuron, i.e. that drift models the mean input to the neuron and diffusion the variance of input?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I agree. You need to provide the pertaining equations here. Are these factors K_ab w_ab and K_ab w_ab^2 from (28) and (29) in the paper, respectively?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for the suggestions. We added a detailed explanation. Note that these factor are quite general and not necessarily in the form of (28) and (29) in the paper.

set_weight( double w )
{
throw BadProperty(
"Please uses the parameters \"drift_factor\" and "

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Typo: "uses" --> "use"

// If the parameter weight is set, we throw a BadProperty
if ( d->known( names::weight ) )
throw BadProperty(
"Please uses the parameters \"drift_factor\" and "

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Typo: "uses" --> "use"

{
throw BadProperty(
"Please uses the parameters \"drift_factor\" and "
"\"diffusion_factor\" to specifiy the weights" );

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Add '.' at end of sentence

if ( d->known( names::weight ) )
throw BadProperty(
"Please uses the parameters \"drift_factor\" and "
"\"diffusion_factor\" to specifiy the weights" );

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Add '.' at end of sentence

non-linearity given by the gain function of the
leaky-integrate-and-fire neuron with delta or exponentially decaying
synapses. The model can be used for a mean-field analysis of spiking
networks.

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Perhaps, a reference to a paper with the actual equations would be useful here?
For instance:
Fourcaud, Nicolas and Brunel, Nicolas, "Dynamics of the firing probability of noisy integrate-and-fire neurons", 14 (2002), pp. 2057--2110.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We added the reference.

double mu,
double sigma )
{
double alpha = 2.0652531522312172;

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I assume this is the evaluation of the Riemann zeta function. Can you add a line comment explaining this briefly?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the suggestion. I added an explanation of the variable alpha


# set up rate neuron and devices
self.rate_neuron = nest.Create(
'lin_rate_ipn', params=self.neuron_params)

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This test only includes the input noise neuron. It would be straightforward to extend this to test the output noise neuron as well. Wouldn't this be a good idea?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, this is a good idea. We extended the test to output noise.


@nest.check_stack
class RateNeuronCommunicationTestCase(unittest.TestCase):
"""Check rate_neuron"""

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You could extend this documentation by giving some more details, such as:
This test checks the correct function of rate neuron with input noise and three transfer functions. It also checks whether delay and weight of the rate connections work correctly.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We updated the documentation according to your suggestions.

The following parameters can be set in the status dictionary.

rate double - Rate (unitless)
tau double - Time constant in ms.

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can you explain the influence of the parameter tau in the model description above? The same for the other transfer functions, tanh and threshold_linear.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I made the documentation a bit more specific stating that it is the time constant of rate dynamics.

@janhahne
Copy link
Contributor Author

@mschmidt87 Thank you for the first review and the valuable suggestions. So far I fixed the typos and the sli unittest. We will have a closer look on the other issues soon.

@heplesser heplesser added this to the NEST 3.0 milestone Jun 29, 2017
Copy link
Contributor

@heplesser heplesser left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@janhahne Impressive work! I have a number of comments including some issues on naming and design that would be good topics for the NEST Open Developer VC (@terhorstd).



/* BeginDocumentation
Name: delay_rate_connection - Synapse type for rate connections with delay.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@janhahne I am not really happy with the model name: It should be delayed_rate_connection. Or, to group things better in the list of models, we should place delayed at the end. For symmetry, we could then have

  • rate_connection_delayed
  • rate_connection_instantaneous
    @terhorstd May you could discuss this in the Open Developer VC?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@heplesser I discussed this with @jschuecker and @ddahmen and we like and apply your second suggestion for renaming.


namespace nest
{

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remove empty line.

* Class representing a delay-rate-connection. A delay-rate-connection
* has the properties weight, delay and receiver port.
*/

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remove empty line.

*/

template < typename targetidentifierT >
class DelayRateConnection : public Connection< targetidentifierT >
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If you decide to change model names, adjust class names accordingly.

class DelayRateConnection : public Connection< targetidentifierT >
{

double weight_;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Remove empty line
  • In NEST code, we generally place data members at the end, even though there are arguments for placing them first. Please move to the end for consistency.
  • Add a doxygen comment for consistency, even though the name is quite clear.

// propagate rate to new time step (exponential integration)
double drive = siegert( P_.tau_m_ * 1e-3,
P_.tau_syn_ * 1e-3,
P_.t_ref_ * 1e-3,
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why do you need to scale the time values here?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good point. The output of the function siegert is a rate in units of Hz. The NEST convention of time constants in ms would therefore give kHz for rates. We improved the conversion and shifted it to the return statement of the siegert functions.

P_.theta_,
P_.V_reset_,
B_.drift_input_[ lag ],
std::sqrt( B_.diffusion_input_[ lag ] ) );
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why the std::sqrt() in the call here? Couldn't you just pass the variance?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We followed your suggestion and pass the variance and convert it to the standard deviation later.


tanh_rate is an implementation of a non-linear rate model with either
input (tanh_rate_ipn) or output noise (tanh_rate_opn) and gain function
Phi(h) = tanh(g * h) and Psi(h) = h for linear_summation = True
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why can the user only choose the gain, but not the location of the inflection point (tanh(gh-x))?.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We followed your suggestion, but implemented the function
tanh(g*(h-x))

@@ -247,6 +250,16 @@ class Event
void set_weight( weight t );

/**
* Set drift_factor of the event (see DiffusionEvent).
*/
virtual void set_drift_factor( weight t ){};
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we use drift_weight and diffusion_weight as names, maybe?�
@terhorstd Discuss in VC.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We decided in the NEST Open Developer VC to keep the old names.

/**
* Event for rate model connections without delay.
*/
class RateNeuronEvent : public SecondaryEvent
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How much of the code here could be moved to SecondaryEvent? There seems quite a bit of overlap with other of these classes.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@heplesser Yes, there is actually a lot of overlap and I am not particular happy about that. I already thought about how to reduce this overlap, but this is not that easy. Each class contains important individual static member variables, which are used in almost every function.

Let me explain in more detail: GapJunctionEvent,RateNeuronEvent,DelayRateNeuronEvent and DiffusionEvent all contain around 200 lines of code, but are more or less identical classes: All of them store arbitrary data of type std::vector<double>. Therefore it would be much better to have one base class (e.g. DoubleVectorEvent) of which all the other classes would be derived.

I could think of two approaches to achieve this:

1. use derived classes
The problem here is that derived classes do not have their own instance of the static variables. One possible solution ( http://stackoverflow.com/a/1390922 ) is to add a virtual function type& get_static() to the base class and call the static variable only with this function. Then each derived class can define its own static variables and by overwriting the virtual get_static function(s) they are also accessible for the functions of the base class.

This approach however does not work for us, as we also have a static function static void set_syn_id( const synindex synid ) in event classes (and in the potential base class) which calls one of the static variables. Calling this static variable with the get_static function yields an cannot call member function without object-error - which means that get_static would have to be static as well, which is not possible as it is a virtual function.

  1. use a template class
    Employ a template <typename something> class DoubleVectorEvent with unused (but unique) template parameters. In this case each template class has its own static members. So one could define the events like typedef DoubleVectorEvent<GapJunction> GapJunctionEvent which would be nice as well. One (however solvable) problem here is the initialization of the static members. A much bigger problem is the operator void DoubleVectorEvent::operator()() which has to stay in event.cpp due to dependencies to node.h. A template definition of the class would however require to move the operator to event.h.

So for now we saw no solution to reduce the overhead here. If you have any suggestion they are of course very welcome!

@janhahne
Copy link
Contributor Author

@mschmidt87 @heplesser Thank you for the careful and detailed review. We finally addressed all of your comments. Please have another look.

@heplesser
Copy link
Contributor

@janhahne Thank you for addressing my comments! I am still not entirely happy with the lack of encapsulation in event handling (see discussion on siegert_neuron()), but I think the code is ready for integration with NEST master. Could you create a follow-up issue to remind us about reviewing the encapsulation issue later?

@jougs Since this is a major addition, I think it makes sense that a non-reviewer merges---would you do the honors?

@janhahne
Copy link
Contributor Author

janhahne commented Aug 9, 2017

@heplesser I created a follow-up issue, see #806

@heplesser @mschmidt87 Thanks again to both of you for the detailed review!

@jougs jougs merged commit 2ffd50c into nest:master Aug 9, 2017
@janhahne janhahne deleted the rate-neuron-framework branch August 9, 2017 17:46
@jougs
Copy link
Contributor

jougs commented Aug 9, 2017

Due to this merge, master is currently failing. Probably it is not not related to the merge itself, but to the length of cppcheck output. See #798 for an explanation and a fix.

@heplesser
Copy link
Contributor

@jougs From the Travis logs, it looks indeed like a tests hitting a resource limit durings static code checking, i.e., the same problem as in #798.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
I: Internal API Changes were introduced in basic internal workings of the simulator that developers need to know S: High Should be handled next T: Enhancement New functionality, model or documentation ZC: Kernel DO NOT USE THIS LABEL ZP: PR Created DO NOT USE THIS LABEL
Projects
None yet
Development

Successfully merging this pull request may close these issues.

6 participants