Build the NwbSortingExtractor spike vector with a single bulk read#4662
Build the NwbSortingExtractor spike vector with a single bulk read#4662h-mayorquin wants to merge 5 commits into
Conversation
h-mayorquin
added
extractors
| # Stable sort by sample_index keeps unit_index ordering on ties, matching the way the generic | ||
| # BaseSorting builder constructs the vector (unit order within a sample). | ||
| order = np.argsort(sample_indices, stable=True) | ||
| n = sample_indices.size |
There was a problem hiding this comment.
Maybe n_spikes to be more explicit?
There was a problem hiding this comment.
Yes, good point.
|
|
||
| # Stable sort by sample_index keeps unit_index ordering on ties, matching the way the generic | ||
| # BaseSorting builder constructs the vector (unit order within a sample). | ||
| order = np.argsort(sample_indices, stable=True) |
There was a problem hiding this comment.
The data isn't presorted in NWB? Or can't be guaranteed to be sorted?
There was a problem hiding this comment.
Yes, we can't guarantee that I added a comment.
| return self._time_vector[frame] | ||
| return frame / self._sampling_frequency + self._t_start | ||
|
|
||
| def _times_to_frames(self, spike_times) -> np.ndarray: |
There was a problem hiding this comment.
in nwb terminology do you say frames instead of samples?
There was a problem hiding this comment.
Yes, nwb favors samples. Do you think it should be samples here? I thought about it and I remember that the get_traces uses frame terminology but maybe we were moving away from it?
There was a problem hiding this comment.
I vote for samples. The spike vector uses spikes["sample_index"] and samples is found across the codebase more often than frames.
There was a problem hiding this comment.
Ok, I changed this.
| """ | ||
| spike_times = np.asarray(self.spike_times_data[:], dtype="float64") | ||
| # boundaries[u]:boundaries[u + 1] is the slice of spike_times belonging to unit u. | ||
| boundaries = np.concatenate(([0], np.asarray(self.spike_times_index_data[:], dtype="int64"))) |
There was a problem hiding this comment.
This seems super fancy and I'm having a hard time following this. So you add a 0 to the beginning then take the diff to get boundaries? Could you give me the dumbed down version of why this works?
There was a problem hiding this comment.
Actually yes, this is kind of complicated. Mmm, let's discuss it tomorrow and we probably should add documentation in hdmf that we c an link for this.
There was a problem hiding this comment.
We discussed this in person and I added a comment.
|
On some real data: https://dandiarchive.org/dandiset/000939/0.260512.1701/files?location=sub-A1808&page=1 |
3 participants
NwbSortingExtractornow overrides_compute_and_cache_spike_vectorto build the whole spike vector from one read of the Units table instead of the generic per-unit loop it inherits fromBaseSorting. An NWB (Neurodata Without Borders) Units table stores every unit's spikes in a single flatspike_timesdataset grouped by unit, with the per-unit boundaries inspike_times_index, so the entire table can be read in one shot and each spike's unit index recovered from those boundaries. The generic builder instead reads one spike train per unit, which is one streamed read per unit; for a remote or streamed file that is the difference between a single request and one request per unit. Construction stays lazy as before: nothing is read until the spike vector is first requested.The seconds-to-samples conversion stays on the segment, shared by the per-unit and bulk paths through
_times_to_frames, and the override only assembles and sorts the structured array, so its output is byte-for-byte identical to the generic builder. Correctness is therefore already exercised by the existing NWB (Neurodata Without Borders) roundtrip tests, and a comment on the override records why it should not be reverted to the generic per-unit path.