Issue with Poisson spikes

I am facing this issue with poisson spiking. The amplitude spikes keep exploding after intervals.


pre = nengo.Ensemble(n_neurons=1, dimensions=1, encoders=[[1]], intercepts=[0], max_rates=[50], neuron_type=nengo.PoissonSpiking(nengo.LIFRate(tau_ref=2/1000)))


p_pre = nengo.Probe(pre.neurons)
plt.plot(sim.trange(),[p_pre]/1000, c=‘k’, label=‘spikes’)

Any idea how do I fix this?

This looks like expected behaviour to me. By ‘exploding’ do you mean those two time-steps where it outputs 2 spikes rather than 1 spike during those time-steps? Assuming you are using a dt of 1ms, given a Poisson distribution with a rate of 50 Hz, this is roughly what one would expect. Specifically, across 2,000 time-steps you would expect ~95 of those steps to contain one spike, ~2 of those steps to contain two spikes, and overall 40x as many steps will contain one spike versus two spikes. (*)

This may be confusing as you are basing the Poisson neuron on the response curve of a LIFRate neuron with a refractory period of 2 ms. LIFRate(tau_ref=0.002) cannot spike more than once in a single millisecond, but when you convert it to Poisson you lose the dynamics of that refractory since a Poisson process does not have a refractory period. Nevertheless, you keep the same response curve, and all this means is that occasionally the neuron spikes more than once in a one millisecond window.

(*) See code checking this mathematically, below:

import scipy.special

def pdf(lmbda, k):
    """Probability density function for a Poisson distribution."""
    return lmbda ** k * np.exp(-lmbda) / scipy.special.factorial(k)

dt = 0.001
rate = 50
one_spike = pdf(rate * dt, k=1)
two_spike = pdf(rate * dt, k=2)

n_steps = 2000
print(one_spike * n_steps)
print(two_spike * n_steps)
print(one_spike / two_spike)


1 Like

Yes, thanks a lot. The two spikes were what I meant by exploding.
I understand the cause now. But cannot think of how to fix this. I really want refractory of 2ms and random spikes so Poisson. Any idea how do I get that? Would decreasing the time step of simulation help ?

Decreasing the time-step would help prevent the case where two spikes are happening in a single step, but there is still always going to be a finite probability of it happening (that probability will just get ‘really’ small for small dt). But either way this wouldn’t give you a refractory period of 2ms. That is, even if you set the time-step to be really small and it spiked at most once per step, you will still be able to find 2ms intervals of time that contain more than one spike.

Generally speaking the two features that you want: (1) a Poisson process, and (2) no two spikes occurring within 2ms of each other… are competing features in a sense. A Poisson process is said to be memoryless in that the system has no ‘memory’ of when it last spiked. In other words, the time of the last event does not affect the probability of the next event. But a refractory period requires a memory of the last event, and therefore such spike trains cannot be described purely as Poisson.

In order to get around this dilemma I think we need to be specific about the features you need. Here are two possibilities depending on what the focus is:

  1. If the refractory period is the important feature, then you could “silence” the process for the 2ms following any spike. It would no longer be Poisson during the 2ms following any given spike, but would give you a refractory period. You could implement this with a custom neuron model, but should be aware that this would alter the response curve from LIFRate by reducing the expected firing rate.
  2. If randomness is the important feature and the type of randomness is not important (i.e., it does not necessarily have to be Poisson), then you could inject random noise into the nengo.LIF neuron model, either by supplying some noise process to the noise parameter of the nengo.Ensemble, or by adding noise to the input signal that you are encoding into the ensemble.

We can help you with whichever path, either way. Two would be a bit easier I think, but both are definitely doable. :slight_smile:

1 Like

Oh just thought of a third option that is related to the second option:

  1. You could use nengo.StochasticSpiking(nengo.LIFRate(tau_ref=2/1000)) with a time-step of dt=0.002. That should give you a randomly varying spike train that does not emit more than one spike per 2ms step.

This option seems the cleanest to me if you don’t need it to be strictly Poisson, but just need some type of randomness. We can go into some details of why this works based on where we go with your replies!

1 Like

Thanks a lot for all sorts of options. Really useful. I wanted to characterise a Hebbian learning rule in nengo. Therefore I wanted all possible spike time difference given fixed frequency and refractory period. The third option sounds easy and sufficient. I shall get back to you after trying this.
Thanks again :smile:

1 Like

Glad I could help. Since you mentioned you are working on a learning rule, another possibility might be to build whatever logic you want regarding the refractory period into the rule itself. For instance, if two spikes occur from the same neuron less than 2ms apart, then you could do something like ignore the later one. Just a thought. Good luck! Let us know how it goes.

Thanks a lot. The solution you suggested earlier is working so I think I can skip adding complexity to learning rule. I am getting good results in characterization. I was thinking of an application to benchmark the unsupervised learning rule similar to BCM, and STDP. I am trying running the examples provided in nengo. Do let me know if you have any suggestions on a better application to benchmark the learning rule.

@arvoelke For the StochasticSpiking solution, why do you need dt=0.002? My understanding is that if the base neuron type never has a rate > 1 / dt, then StochasticSpiking will never emit more than one spike per timestep. For any LIFRate with tau_ref >= 1 / dt, you’ll never get more than one spike per timestep, so for tau_ref = 2 / 1000, any dt <= 0.002 will work (i.e. the default dt=0.001 would also work fine).

Because I was addressing this feature:

That is, to satisfy the refractory period constraint. Not whether you get more than one spike per timestep. Rather, whether two spikes occur within 2ms of each other.

You are right that if we only care about the multiple spikes per step, then any dt <= 0.002 will work to prevent that with nengo.StochasticSpiking(nengo.LIFRate(tau_ref=0.002)). But if you additionally want to ensure that they are 2ms apart then you need the 2ms step. For instance with dt = 0.001 you could get two spikes that are 1ms apart, and with dt = 0.0001 you get two spikes that are 0.1ms apart (albeit with very low probability).

That makes sense. When I had read “I really want refractory of 2ms” I was just thinking of it with respect to the tuning curve (i.e. saturating at 500 Hz), not with respect to the spike train.

One way to do spike trains that don’t allow spikes on adjacent timesteps would be to generate the inter-spike intervals (ISIs) randomly, and then form the spike train based on the ISIs. We don’t currently have a neuron type that does this, but it wouldn’t be too hard to implement (our current framework should support it).

1 Like

Thanks. I actually wanted the ISI to be greater than 2ms. I am not sure if increasing the time step would be good, because I want this refractory for only one ensemble and not all.