SpikeGLXRawIO incorrectly handles folder per probe organization.

[grahamfindlay]

Describe the bug
SpikeGLX may write data using "folder per probe" organization. Here is the description, from the SpikeGLX User Manual:

Note that there is an option on the Save Tab called Folder per probe.
If this is set, there is still a run folder for each g-index run-name_gN.
However, inside that there is also a subfolder for each probe that contains
all the t-indices for that g-index and that probe. A probe subfolder is
named like run-name_gN/run-name_gN_imecM.

For example:

- run_dir/ (example: 3-1-2021/)
  - gate_dir_g0/ (example: 3-1-2021_g0/)
  ...
  - gate_dir_gN/
    - probe_dir_imec0/ (example: 3-1-2021_g0_imec0/)
    ...
    - probe_dir_imecN/
      - trigger_file_gN_imecN_t0.lf.bin
      - trigger_file_gN_imecN_t0.lf.meta
      - trigger_file_gN_imecN_t0.ap.bin
      - trigger_file_gN_imecN_t0.ap.meta
      ...
      - trigger_file_gN_imecN_tN.lf.bin
      - trigger_file_gN_imecN_tN.lf.meta
      - trigger_file_gN_imecN_tN.ap.bin
      - trigger_file_gN_imecN_tN.ap.meta

This issue applies when using the folder-per-probe layout and reading probes with a device index greater than 0. For example, trying to read the following file:

3-1-2021/3-1-2021_g0/3-1-2021_g0_imec1/3-1-2021_g0_imec1.ap.bin

will result in the following error:

AssertionError: stream_id imec1.ap is not in [np.str_('imec0.ap'), np.str_('imec0.lf')]

It should be checking for membership in [np.str_('imec1.ap'), np.str_('imec1.lf')], as that is the only probe in the directory.

The incorrect stream IDs parsed for this directory are generated in neo.rawio.spikeglxrawio.scan_files(), on line 405.

stream_name = f"{device_kind}{device_index}{stream_kind}"

This creates the string {imec}{0}{.ap} because device_index is 0, when it should be 1.

The reason that device_index is 0 is because of lines 388-398 in neo.rawio.spikeglxrawio.scan_files():

    # TODO: handle one box case
    info_list_imec = [info for info in info_list if info.get("device") != "nidq"]
    unique_probe_tuples = {get_probe_tuple(info) for info in info_list_imec}
    sorted_probe_keys = sorted(unique_probe_tuples)
    probe_tuple_to_probe_index = {key: idx for idx, key in enumerate(sorted_probe_keys)}

    for info in info_list:
        if info.get("device") == "nidq":
            info["device_index"] = ""  # TODO: Handle multi nidq case, maybe use meta["typeNiEnabled"]
        else:
            info["device_index"] = probe_tuple_to_probe_index[get_probe_tuple(info)]

The problem is that you cannot enumerate the sorted probe keys on line 392 and expect that enumeration index to be the correct device index on line 398. In the case of folder per probe organization, where the the probe is imecN, the enumeration index will always be 0, but the device index should be N.

One way to solve this could be to parse the device index directly from the probe name or the info["device"] field, both of which should be "imecN":

info["device_index"] = int(re.split(r"(^[^\d]+)", info["device"])[-1])

If this index cannot be assigned directly to info["device_index"], it could at least be saved and used to offset the probe index from the probe_tuple_to_index dict. Someone with access to a wider variety of SpikeGLX outputs (non-folder-per-probe layouts, OneBox, etc) than I do should double-check!

This bug appears to have been introduced in #1608.

Environment:

• OS: Linux
• Python version 3.13.2
• Neo version 0.14.0
• NumPy version 2.1.3

Additional context
One other thing I noticed while tracking down this bug: On lines 530-534 in neo.rawio.spikeglxrawio.extract_stream_info, device and stream_kind are parsed twice each in quick succession: first from a regex, and then by just splitting the filename. Since they're obtained two different ways, might not hurt to assert that they match. Or just not bother recalculating them the second time, since they're already in hand? Nbd, just something that caught my eye.

Thank you!

[zm711]

@h-mayorquin and @alejoe91, you are probably best to handle this.

[h-mayorquin]

Hi, @grahamfindlay

Would it be possible for you to share:

1. The exact call that triggers the error (are you calling the top level folder or something else?)
2. Can you share example data with us? Maybe you can stub it with something like this:
   https://gist.github.com/h-mayorquin/b5abe6ae8e8484976c2d27e0bb04e2a2

Or can it be reproduced with any of the datasets that you previously shared?
https://gin.g-node.org/NeuralEnsemble/ephy_testing_data/src/master/spikeglx/multi_trigger_multi_gate

That way, we would avoid a regression like this again.

[grahamfindlay]

Hi, @h-mayorquin

The error can be reproduced with several of the datasets I previously shared. Specifically, any of the ones that I labeled as using folder-per-probe organization in the README (i.e., -CI0, -CI2, and -CI4). See below for a minimal reproducible example that uses -CI0.

The exact call that triggers the error is actually in spikeinterface. Sorry, I totally failed to mention that. Example:

from pathlib import Path
import spikeinterface.extractors as se

run_dir = Path("./ephy_testing_data/spikeglx/multi_trigger_multi_gate/SpikeGLX/5-19-2022-CI0")
gate_dir = run_dir / "5-19-2022-CI0_g0"
probe_dir = gate_dir / "5-19-2022-CI0_g0_imec1"
assert probe_dir.exists()

extractor = se.SpikeGLXRecordingExtractor(gate_dir, stream_id="imec1.ap")

The error is thrown in `spikeinterface.extractors.neoextractors.NeoBaseRecordingExtractor, line 218:

assert stream_id in stream_ids, f"stream_id {stream_id} is not in {stream_ids}"

The reason is that stream_ids is taken from self.neo_reader.header["signal_streams"]["name"], and is array(['imec0.ap', 'imec0.lf'], dtype='<U64').

So, basically, neo will happily create a reader and parse the header (incorrectly) without throwing an error, but then spikeinterface raises the error when it tries to validate the user arguments against the incorrectly parsed metadata.

reader = neo.rawio.SpikeGLXRawIO(gate_dir)
reader.parse_header() # No error thrown, but this is where the root of the problem actually lies. 

So I don't think any changes should be required downstream in spikeinterface.

Hope that helps!
-Graham

[zm711]

Okay, so in this case I think we will also need to add a new test into our ci that parses the possible streams based on the files present and then ensures they actually exist. This is good to know Graham thanks.

[h-mayorquin]

Hi, @grahamfindlay

This works just fine for me:

from pathlib import Path

folder_path = Path("/home/heberto/neuroconv_testing_data/ephy_testing_data/spikeglx")
assert folder_path.exists()

folder_path  = folder_path / "multi_trigger_multi_gate/SpikeGLX/5-19-2022-CI0"
assert folder_path.exists()

gate_dir = folder_path / "5-19-2022-CI0_g0"
probe_dir = gate_dir / "5-19-2022-CI0_g0_imec1"

from spikeinterface.extractors import SpikeGLXRecordingExtractor

recording = SpikeGLXRecordingExtractor(gate_dir, stream_id="imec1.ap")
recording

And it seems it detects the streams it should:

Can you help me see if I am making am mistake.

[grahamfindlay]

This had me very confused for a minute!

It apparently matters whether the other probe is present: to load imec1, imec0 also has to be present. If you just gin get-content 5-19-2022-CI0_g0_imec1, you will run into the error. But if you then do gin get-content 5-19-2022-CI0_g0_imec0, it will go away.

[grahamfindlay]

Incidentally, it also appears that you can elicit this error even if you have both probes downloaded and present in sibling directories, by passing the probe_dir instead of the gate_dir to the SpikeGLXRecordingExtractor, like so:

from pathlib import Path
import spikeinterface.extractors as se

run_dir = Path("./ephy_testing_data/spikeglx/multi_trigger_multi_gate/SpikeGLX/5-19-2022-CI0")
gate_dir = run_dir / "5-19-2022-CI0_g0"
probe_dir = gate_dir / "5-19-2022-CI0_g0_imec1"
assert probe_dir.exists()

extractor = se.SpikeGLXRecordingExtractor(probe_dir, stream_id="imec1.ap")

From the documentation of the SpikeGLXRecordingExtractor class, it seems like this should be allowed:

    folder_path : str
        The folder path to load the recordings from.
    ...
    stream_id : str or None, default: None
        If there are several streams, specify the stream id you want to load.
        For example, "imec0.ap", "nidq", or "imec0.lf"
    ...
    Examples
    --------
    >>> from spikeinterface.extractors import read_spikeglx
    >>> recording = read_spikeglx(folder_path=r'path_to_folder_with_data', load_sync_channel=False)

None of which implies that you couldn't just pass the probe folder to the constructor.

It's sort of the same situation with the docs for SpikeGLXRawIO:

    Parameters
    ----------
    dirname: str, default: ''
        The spikeglx folder containing meta/bin files

Which would actually cause me to try passing the probe_dir before trying the gate_dir, although shortly afterwards the same docstring does mention subfolders, which implies the gate_dir:

    * This IO reads the entire folder and subfolders locating the `.bin` and `.meta` files

[h-mayorquin]

OK, I managed to reproduce this by moving only the folder containing the data from one probe and one gate and isolating it from the rest of your data:

from pathlib import Path

folder_path = Path("/home/heberto/data/spikeglx_only_one_probe")
assert folder_path.exists()

probe_folder_path = folder_path / "5-19-2022-CI0_g0_imec1"

from spikeinterface.extractors import read_spikeglx

read_spikeglx.get_streams(probe_folder_path)

recording = read_spikeglx(probe_folder_path, stream_id="imec1.ap")
recording

and as you mention if you pass the data containing only one probe that triggers the error as well:

from pathlib import Path

folder_path = Path("/home/heberto/neuroconv_testing_data/ephy_testing_data/spikeglx")
assert folder_path.exists()

folder_path  = folder_path / "multi_trigger_multi_gate/SpikeGLX/5-19-2022-CI0"
assert folder_path.exists()

gate_dir = folder_path / "5-19-2022-CI0_g0"
probe_dir = gate_dir / "5-19-2022-CI0_g0_imec1"

from spikeinterface.extractors import SpikeGLXRecordingExtractor

recording = SpikeGLXRecordingExtractor(probe_dir, stream_id="imec1.ap")
recording

So with that we can work, but now that I have understood better I do have the question:

Is it that common that the data for probe 1 would be obtained without the accompanying data for probe 0? It would be good for me to understand how common this is before moving forward.

[h-mayorquin]

The problem is that you cannot enumerate the sorted probe keys on line 392 and expect that enumeration index to be the correct device index on line 398. In the case of folder per probe organization, where the the probe is imecN, the enumeration index will always be 0, but the device index should be N.

One way to solve this could be to parse the device index directly from the probe name or the info["device"] field, both of which should be "imecN":

This is the way it was before, but the problem is that sometimes people change the names of their files. Check this out:

#1606

There is an entry in the meta called fileName that keeps the original name. Maybe it would be as simple as using that one.

The fundamental problem is that the meta file does not contain the device_index. The calculation I'm using is here:

python-neo/neo/rawio/spikeglxrawio.py

Lines 379 to 386 in f05b92f

	# Probe index calculation
	# The calculation is ordered by slot, port, dock in that order, this is the number that appears in the filename
	# after imec when using native names (e.g. imec0, imec1, etc.)
	def get_probe_tuple(info):
	slot = normalize(info.get("probe_slot"))
	port = normalize(info.get("probe_port"))
	dock = normalize(info.get("probe_dock"))
	return (slot, port, dock)

This approach is from @jenniferColonell, but it breaks down when you don't know if there were other ports, docks, or slots available, since those numbers aren't available when looking at the meta of a single probe.

I think the best solution is to add the fix using the original filename, but I'm a bit concerned that we're exchanging one type of brittleness (parsing device_index from the filename) for another (calculating it from all the available metadata). Both approaches have their failure modes, but I think we broke backward compatibility when we solved #1606, and it would be good to undo that regression.

I will briefly discuss this in the meeting tomorrow and then push a fix. I've already got a test working on the CI with our current data. Let me know what do you think @grahamfindlay.

[jenniferColonell]

I'm pretty sure I don't use slot/port/dock values anywhere, because they aren't generally important in the output data! Although they are included in the metadata (to trace hardware issues) they aren't meant to tell you the probe index for the run. Parsing the original filename is likely the best bet for handling this case -- unless it is actually data from the 3A era.
It's true -- if the user has changed the file names AND the name in the metadata, there will be an issue. But reading the data is vulnerable to all kinds of uncontrolled changes to the metadata.

[grahamfindlay]

I think that the first case, where folder-per-probe organization was/is used to create the directory structure/names, but the data from the two probes actually do not exist in sibling directories, is probably very uncommon. In my case, it happened because one probe malfunctioned but continued to acquire data, and rather than keep the many terabytes of junk around, I just deleted the bad probe's data.

Being able to pass the probe_dir instead of the gate_dir on the other hand, probably does make sense to support (or documentation could be changed to make it explicit that gate_dir is required, but this may be confusing or surprising for some users). Supporting probe_dir could be as easy as just detecting it using Path.parent to get the gate_dir. However, this would still maintain the limitation that you can't really load one probe without having the other around. So I think, ideally, it would be possible to just load one probe without parsing its sibling's meta files. After all, there isn't a fundamental reason I can think of why you should really have to touch, peek at, or even know about imec0 to load imec1 and work with it.

I see what you mean about device index not being present in the metadata apart from the filename, and about trading one kind of brittleness for another. If it were me, and I had to pick, I would decline to support arbitrary file renaming, including stripping the device index, unless there was a really good reason for users to be doing it, especially if they could be un-doing it. I suppose you could check if the device index is present in the filename, and if not, fall back on the current method. I think the best solution sounds like the one you proposed: using the fileName meta entry. If someone's going and arbitrarily renaming that, they're a psychopath! I don't think the official SpikeGLX suite of tools will even tolerate files whose names don't match what's in their metadata -- it's considered a failure of validation. But don't quote me on that.

[h-mayorquin]

I'm pretty sure I don't use slot/port/dock values anywhere, because they aren't generally important in the output data! Although they are included in the metadata (to trace hardware issues) they aren't meant to tell you the probe index for the run. Parsing the original filename is likely the best bet for handling this case -- unless it is actually data from the 3A era. It's true -- if the user has changed the file names AND the name in the metadata, there will be an issue. But reading the data is vulnerable to all kinds of uncontrolled changes to the metadata.

Thanks for jumping in.

What I meant is that I once asked how to calculate the probe_index that is included in the filename here

Aug 2nd, 2024 at 9:34 AM
The final number is the probe index for the run; the probe indices are ordered according to slot,port,and dock (for NP 2.0 probes). The probe index is also how the probes are referenced in the SGLX gui.

Just to make that more explicit: imec0 is the probe in the lowest slot and port number, imec1 in the next highest, and so on. If you have probes in {slot 2, port 3}, {slot 3, port1} and {slot3, port2}, the probe indices in the SGLX output will be 0, 1, and 2, respectively.

Probably my misunderstanding somewhere there.

[h-mayorquin]

I think that the first case, where folder-per-probe organization was/is used to create the directory structure/names, but the data from the two probes actually do not exist in sibling directories, is probably very uncommon. In my case, it happened because one probe malfunctioned but continued to acquire data, and rather than keep the many terabytes of junk around, I just deleted the bad probe's data.

Being able to pass the probe_dir instead of the gate_dir on the other hand, probably does make sense to support (or documentation could be changed to make it explicit that gate_dir is required, but this may be confusing or surprising for some users). Supporting probe_dir could be as easy as just detecting it using Path.parent to get the gate_dir. However, this would still maintain the limitation that you can't really load one probe without having the other around. So I think, ideally, it would be possible to just load one probe without parsing its sibling's meta files. After all, there isn't a fundamental reason I can think of why you should really have to touch, peek at, or even know about imec0 to load imec1 and work with it.

Thanks, @grahamfindlay that sounds reasonable to me. I just opened #1661 which should fix the error here.

Thanks for working with me through this fix.