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

Fix postsynaptic trace edge case #1137

Merged
merged 47 commits into from
Apr 8, 2019
Merged

Conversation

clinssen
Copy link
Contributor

This fixes #1034.

With props to my fellow sleuth @suku!

Background:

The synapse, when handling an incoming presynaptic spike at time t_pre_spike in send(), requests postsynaptic history within the interval t_prev_pre_spike - dendritic_delay < t ≤ t_pre_spike - dendritic_delay. Spikes in this interval will have their access counters incremented. The postsynaptic neuron, in Archiving_Node::set_spiketime(), will purge all spikes from the buffer which have been accessed by all connecting synapses.

A minimal reproducing example of the bug has presynaptic spike times at [2, 3, 10] and postsynaptic at [1, 2, 3] (all times in ms). All delays (axonal, dendritic, synaptic) are equal to 1. The failure occurs as follows: The first postsynaptic spike at t = 2 is added to the postsynaptic buffer. The first presynaptic spike at t = 3 causes the synapse to access the postsynaptic buffer. This increments the access counter of the spike at t = 2, marking it for deletion. An extra postsynaptic spike at t = 3 is necessary because the history always needs to have at least two entries if a spike is to be deleted. Next, the postsynaptic neuron spikes at t = 4, which causes the spike to be deleted from the buffer. When another presynaptic spike then arrives at t = 4, the get_K_value() function is called on the postsynaptic partner to obtain the trace at t = 3, but the necessary spike at time 2 has now been deleted. Consequently, zero is erroneously returned for the trace value.

How this PR fixes things:

Recall that spike buffer pruning occurs upon firing of the neuron at time t_sp in Archiving_Node::set_spiketime(). A new condition for pruning is added: only remove a spike if there is another, later spike, that is strictly more than delay away from t_sp (https://github.com/clinssen/nest-simulator/blob/76a71602edc929116a464d4e8ab34c8cc00b5242/nestkernel/archiving_node.cpp#L210). This ensures that when get_K_value() is called to obtain the trace at time t - delay, a spike strictly earlier than this is still in the buffer.

To obtain the maximum delay value across all synapses, a new private member variable max_delay_ is added in Archiving_Node, which is updated whenever Archiving_Node::register_stdp_connection() is called.

Notes for reviewers:

@terhorstd terhorstd added T: Bug Wrong statements in the code or documentation ZC: Infrastructure DO NOT USE THIS LABEL S: High Should be handled next I: No breaking change Previously written code will work as before, no one should note anything changing (aside the fix) ZP: PR Created DO NOT USE THIS LABEL labels Mar 8, 2019
Copy link
Contributor

@hakonsbm hakonsbm left a comment

Choose a reason for hiding this comment

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

Thank you for finding and fixing this bug! I tried benchmarking with the brunel_ps.sli example, and it shows no significant performance loss from caching in trace_. In the C++ part there is some formatting that needs fixing, see below.

The added testfile, test_post_trace.py, is quite messy. Some suggestions:

  • The plotting isn't required for running the unittest, so I suggest to have it separate from the test. perhaps commented out below together with something like a "make_plot" function. It also needs some refactoring, but as it's only for debugging that's not a priority.
  • The helper functions run_post_trace_test_nest_, run_post_trace_test_python_reference_, and nest_trace_matches_ref_trace can be moved out of the PostTraceTestCase class. One possibility is to create a separate class containing them as methods, which can reduce the number of arguments.

, last_spike_( -1.0 )
, Ca_t_( 0.0 )
, Ca_minus_( 0.0 )
, tau_Ca_( 10000.0 )
, beta_Ca_( 0.001 )
, trace_(0.)
Copy link
Contributor

Choose a reason for hiding this comment

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

Variables should be initialized in the same order as they are declared in the class. So trace_ should be initialized directly after max_delay_.

, last_spike_( n.last_spike_ )
, Ca_t_( n.Ca_t_ )
, Ca_minus_( n.Ca_minus_ )
, tau_Ca_( n.tau_Ca_ )
, beta_Ca_( n.beta_Ca_ )
, trace_(n.trace_)
Copy link
Contributor

Choose a reason for hiding this comment

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

Same as above, move initialization of trace_ to directly after max_delay_.

}

// search for the latest post spike in the history buffer that came strictly before `t`
Copy link
Contributor

Choose a reason for hiding this comment

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

Line is too long (>80).

}
return 0;

// this case occurs when the trace was requested at a time precisely at or before the first spike in the history
Copy link
Contributor

Choose a reason for hiding this comment

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

Line is too long (>80).

// except the penultimate one. we might still need it.
// only remove a spike if:
// - its access counter indicates it has been read out by all connected STDP synapses, and
// - there is another, later spike, that is strictly more than (max_delay_ + eps) away from the new spike (at t_sp_ms)
Copy link
Contributor

Choose a reason for hiding this comment

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

Keep line lengths under 80 characters.

#
# 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 thed
Copy link
Contributor

Choose a reason for hiding this comment

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

"thed" -> "the"

import scipy as sp
import scipy.stats
import matplotlib.pyplot as plt
import matplotlib.ticker as plticker
Copy link
Contributor

Choose a reason for hiding this comment

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

math is unused, and sp.stats is used throughout, so no need to import scipy.stats.


nest.set_verbosity("M_WARNING")

post_weights = {'parrot': [], 'parrot_ps': []}
Copy link
Contributor

Choose a reason for hiding this comment

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

post_weights isn't used anywhere.

'delay': delay})

# get STDP synapse
syn_ps = nest.GetConnections(source=pre_parrot_ps,
Copy link
Contributor

Choose a reason for hiding this comment

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

syn_ps isn't used after assigning it.

t = trace_nest_t[i]
print("* Finding ref for NEST timepoint t = "
+ str(t) + ", trace = " + str(trace_nest[i]))
for i_search, t_search in enumerate(
Copy link
Contributor

Choose a reason for hiding this comment

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

Doesn't look like i_search is used anywhere, so no need to enumerate.

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.

@clinssen Please minimize the test to a pure test script with no plotting. The test should not contain anything that is not required for testing. Also, move it to the regressiontest directory and rename it issue-1034. Our test system handles python scripts in all test directories. If you want to keep the plotting for documentation purposes, a notebook under doc/model_details would be a good place.

@heplesser heplesser removed the request for review from jarsi March 25, 2019 12:35
@clinssen
Copy link
Contributor Author

Gentlemen, many thanks for your excellent suggestions. The new series of commits should address all of the points you raised. Additionally, I've been able to tighten the (absolute and relative) tolerances from 1E-2 to 1E-3 by making some modifications in how the reference timeseries is generated.

Copy link
Contributor

@hakonsbm hakonsbm left a comment

Choose a reason for hiding this comment

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

Thanks for the changes, looks a lot better now! I have a couple remarks:

  • I'm not sure if we need a copyright notice in notebooks (there aren't any in the others)
  • It doesn't look like the added NEST name pre_trace is used for anything. It can probably be removed.

I also have some suggestions about the python test, see below.

if ( history_.front().access_counter_ >= n_incoming_ )
const double next_t_sp = history_[ 1 ].t_;
if ( history_.front().access_counter_ >= n_incoming_
&& t_sp_ms - next_t_sp > max_delay_
Copy link
Contributor

Choose a reason for hiding this comment

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

Use and instead of &&.

nest.Connect(
pre_parrot_ps + post_parrot_ps,
spikes
)
Copy link
Contributor

Choose a reason for hiding this comment

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

Why not just nest.Connect(pre_parrot_ps + post_parrot_ps, spikes) in one line?

for sp_idx in range(n_spikes):
t_sp = self.post_spike_times_[sp_idx] \
+ self.delay_ \
+ self.dendritic_delay_
Copy link
Contributor

Choose a reason for hiding this comment

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

It is recommended to use parentheses to break lines instead of backslash.

/ self.tau_minus_)

n_spikes = len(self.pre_spike_times_)
for sp_idx in range(n_spikes):
Copy link
Contributor

Choose a reason for hiding this comment

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

Loop over spike times, not indices.

"""

n_timepoints = len(trace_nest_t)
for i in range(n_timepoints)[1:]:
Copy link
Contributor

Choose a reason for hiding this comment

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

Loop over trace_nest_t, not indices.

break

if (not traces_match) \
and i_search == len(self.pre_spike_times_) - 1:
Copy link
Contributor

Choose a reason for hiding this comment

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

Use parentheses to break lines.

post_spike_times2,
pre_spike_times2]

for spike_times_idx in range(len(pre_spike_times)):
Copy link
Contributor

Choose a reason for hiding this comment

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

Instead of iterating the index, iterate the spike times with zip. So it would become something like

for pre_spike_time, post_spike_time in zip(pre_spike_times,
                                           post_spike_times):

+ ", ".join([str(t) for t in post_spike_times]) + "]")

for delay in delays:
dendritic_delay = delay
Copy link
Contributor

Choose a reason for hiding this comment

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

dendritic_delay is unused.

tau_minus=tau_minus,
trace_match_atol=1E-3,
trace_match_rtol=1E-3)
assert test.nest_trace_matches_python_trace()
Copy link
Contributor

Choose a reason for hiding this comment

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

Use self.assertTrue().



if __name__ == "__main__":
runner = unittest.TextTestRunner(verbosity=99)
Copy link
Contributor

Choose a reason for hiding this comment

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

Verbosity only goes up to 2.

@heplesser
Copy link
Contributor

@clinssen

  • I'm not sure if we need a copyright notice in notebooks (there aren't any in the others)

I would leave the copyright notice out of the notebook for now, in accordance with the existing notebooks.

@heplesser heplesser added this to the NEST-2.18 milestone Apr 1, 2019
@hakonsbm
Copy link
Contributor

hakonsbm commented Apr 5, 2019

@clinssen Thanks for the changes! There are a couple comments in the python test you may have overlooked:

  • L131, use parentheses to break lines here
  • Can the if-tests on lines L192 and L207 be merged together?

I'll give the thumbs up when these comments are resolved.

@heplesser
Copy link
Contributor

@clinssen Thank you for the notebook. Could you add some more text in the final part of the notebook explaining the figures, especially the synaptic trace from NEST? Is the behaviour shown in the graphs the behaviour before or after this fix has been applied? If it is the behaviour after the fix, then why do we still see this difference?

@clinssen
Copy link
Contributor Author

clinssen commented Apr 8, 2019

The latest commits should address all of your comments. @hakonsbm: the if-tests cannot be merged, because I'm using the pattern:

valid = test_if_valid1()
if not valid:
    valid = test_if_valid2()
if not valid:
    ...

@clinssen
Copy link
Contributor Author

clinssen commented Apr 8, 2019

@heplesser: I posted my previous comment on a stale version of this page; please hang tight!

Copy link
Contributor

@hakonsbm hakonsbm left a comment

Choose a reason for hiding this comment

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

@clinssen Fair enough, I won't nitpick any more. 👍

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
I: No breaking change Previously written code will work as before, no one should note anything changing (aside the fix) S: High Should be handled next T: Bug Wrong statements in the code or documentation ZC: Infrastructure 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.

get_K_value() can return 0 erroneously in edge case
5 participants