[Bug] IndexError in spike_train_synchrony.py, annotate_synchrofacts

[Moritz-Alexander-Kern]

This Bug was originally discovered by @skrausse , thank you for reporting.

Describe the bug When trying to use the delete_synchrofacts method from Synchrotool class in spike_train_synchrony.py, in specific cases an IndexError is raised.

To Reproduce

1. Install elephant and neo
2. Run the following minimal example:

# fix synchrotool
import elephant, neo
import quantities as pq
from elephant import spike_train_synchrony


import elephant.spike_train_generation as stg
import numpy as np

np.random.seed(1)  # set random seed 2 for no error

SPP = stg.StationaryPoissonProcess(50 * pq.Hz, t_start=0 * pq.ms,
                                   t_stop=1000*pq.ms)
test_sts = SPP.generate_n_spiketrains(2)

sampling_rate = 30000*pq.Hz
test_obj = spike_train_synchrony.Synchrotool(
    test_sts, sampling_rate=sampling_rate, spread=2)
test_obj.delete_synchrofacts(threshold=2, in_place=True)

Gives the following error messages: (if np.random.seed(1), if the random seed is set to 2 , no error is raised)

Traceback (most recent call last):
    test_obj.delete_synchrofacts(threshold=2, in_place=True)
  elephant/spike_train_synchrony.py", line 325, in delete_synchrofacts
    self.annotate_synchrofacts()
  elephant/spike_train_synchrony.py", line 407, in annotate_synchrofacts
    complexity_per_spike = epoch_complexities[spike_to_epoch_idx]
IndexError: index 91 is out of bounds for axis 0 with size 91

Expected behavior Index should not be out of bounds

Environment

• OS (e.g., Linux):Ubuntu 20.04.4 LTS (64-bit)
• How you installed elephant (conda, pip, source): pip, source
• Python version: 3.8
• neo python package version: 0.10.2
• elephant python package version: 0.11.0
• quantities: 0.13.0
• numpy: 1.22.3

[Moritz-Alexander-Kern]

Hot-fix: Added an if-clause, which checks if index is out of bounds, see commit: 786194c2378ad866776f74c6c510512d39177fd7 https://github.com/INM-6/elephant/tree/fix/synchrotool

In order to test if this Bug still occurs in connection with specific .nix files, the following beta release was created for further testing:

test.pypi: https://test.pypi.org/project/elephant-test/0.11.2b1/

Installation: pip install -i https://test.pypi.org/simple/ elephant-test==0.11.2b1

Update 20.10.2023: the following is no longer up to date, see latest comments ~~This bug seems to be related to rounding in quantities and could be traced down to:~~ https://github.com/NeuralEnsemble/elephant/blob/0df45812b93229b3186c061cb84ec29fd380e09c/elephant/spike_train_synchrony.py#L399-L407

~~More specifically line 406, in certain cases the maximum value here can be larger than the maximum value in right_edges. This could be due to rounding errors.~~

[Kleinjohann]

Hi! In this minimal example, the problem is actually that the spiketrains do not have a sampling rate. The Synchrotool assumes that the spiketrains are actually discrete and that all spike times can only be multiples of the sampling rate that is passed as a parameter. We have a utility function for discretising such artificial spiketrains, when we apply that no error is raised:

# fix synchrotool
import elephant, neo
import quantities as pq
from elephant import spike_train_synchrony
from elephant.conversion import discretise_spiketimes


import elephant.spike_train_generation as stg
import numpy as np

np.random.seed(1)  # set random seed 2 for no error

SPP = stg.StationaryPoissonProcess(50 * pq.Hz, t_start=0 * pq.ms,
                                   t_stop=1000*pq.ms)
test_sts = SPP.generate_n_spiketrains(2)

sampling_rate = 30000*pq.Hz
test_sts = discretise_spiketimes(test_sts, sampling_rate)

test_obj = spike_train_synchrony.Synchrotool(
    test_sts, sampling_rate=sampling_rate, spread=2)
test_obj.delete_synchrofacts(threshold=2, in_place=True)

[Kleinjohann]

This error can however also be caused by

1. (t_stop - t_start) % sampling_rate != 0 which means that not the entire period of the spike train will be binned since otherwise the last bin would be incomplete
2. st[-1] == st.t_stop which would cause the last spike to be missed

I see multiple options for addressing this:

• add input checks asserting that (t_stop - t_start) % sampling_rate == 0 and that st[-1] < st.t_stop
• catch the error and raise a more meaningful one
• keep the problematic spikes out of the calculation, then array_annotate them with -1 or some other error value and raise a warning

Either way we also need to make this issue more clear in the documentation of the function.

[Moritz-Alexander-Kern]

The issue could be traced to the following:

Thanks @skrausse who posted on INM-6/elephant:

[...] @jo460464 already found the relevant code in the complexity class, where the t_stop bin is excluded and therefore leads to a failure.

https://github.com/NeuralEnsemble/elephant/blob/0354d0ed384592f2bf9aee21a559a4a7b4b59d48/elephant/statistics.py#L1418-L1446

Shifting the bins by half a sampling period leaves events at t_stop to be out of bounds.

[skrausse]

Discussion at elephant meeting (23.10.23):

• binning issue traced back synchrotool -> complexity -> time_histogram -> binned_spiketrains
• if provided with sampling rate / period: Shifting of spiketrains and adding an epoch would solve the problem
  • Issue: binned spiketrain relies on shared t_start and t_stop between spiketrains and provided values to not make assumptions on sampling rate
• proposed fix: condition on binned spiketrain: if sampling rate provided and shift parameter allowed, then shift spiketrains; else warn for spikes being removed if that coincides with t_stop

[Moritz-Alexander-Kern]

We also discussed the suggested changes on branch: https://github.com/jo460464/elephant/tree/enh/iss101