[Bug]: EventDetection `source_idx` shape doesn't validate
What happened?
Hello!
I'm attempting to convert my data to the NWB format, and have some voltage traces recorded from a multi-electrode array. Spikes were detected in these traces. If I want to save the traces and the spike events, the way to go seems like I should use the EventDetection data interface.
If my raw data is two-dimensional, e.g., 100,000 samples and 1000 channels, would the source_idx also be two dimensional? The documentation seems to suggest so:
Indices (zero-based) into source ElectricalSeries::data array corresponding to time of event. Module description should define what is meant by time of event (e.g., .25msec before action potential peak, zero-crossing time, etc). The index points to each event from the raw data
However, when I attempt to validate my file with pynwb-validate, here is the output:
Validating div21-carbachol30um.nwb against cached namespace information using namespace 'core'.
- found the following errors:
EventDetection/source_idx (processing/ecephys/mxw_events/source_idx): incorrect shape - expected '[None]', got '(2469376, 2)'
The two columns of the source_idx array are the sample and channel indices into my ElectricalSeries, which has shape (19200200, 993).
Steps to Reproduce
from datetime import datetime
from uuid import uuid4
import numpy as np
from dateutil.tz import tzlocal
from pynwb import NWBHDF5IO, NWBFile
from pynwb.ecephys import ElectricalSeries, SpikeEventSeries
import h5py
from pathlib import Path
datapath = Path("/Data/pharma-spkdet-recs/expts/carbachol/c0087/")
f = h5py.File(datapath / "div21-carbachol30um.raw.h5")
model_name = f['wellplate/version'][0].decode('UTF-8')
serial_number = f['wellplate/id'][0].decode('UTF-8')
mapping = f['data_store/data0000/settings/mapping']
starttime = datetime.fromtimestamp(f['data_store/data0000/start_time'][0]/1000,tz=tzlocal())
try:
platingdate = f['wellplate/well000/Plating Date'][0].decode('UTF-8')
except KeyError:
platingdate = ""
nwbfile = NWBFile(
session_description="30 uM Carbachol added to rat cortical culture",
identifier=str(uuid4()),
session_start_time=starttime,
experiment_description="".join(("Applied 30 uM Carbachol via a pipette to the well of an in vitro rat ",
"cortical culture. Application was at 60s in the recording (no exact timestamp). ",
f"Culture was plated on {platingdate}")),
keywords=["ecephys", "carbachol"]
)
device = nwbfile.create_device(
name="hdmea",
description="1024 Channel High-Density Microelectrode Array",
manufacturer="MaxWell Biosystems",
model_name=f'{model_name}',
serial_number=f'{serial_number}',
)
el_group = nwbfile.create_electrode_group(
name="well0",
description="Electrodes in the single-well of the MaxOne",
device=device,
location="in vitro"
)
for ch,el,x,y in mapping:
nwbfile.add_electrode(
group=el_group,
#label=f"elec{el}",
location="well0",
id=ch,
x=x,
y=y
)
all_table_region = nwbfile.create_electrode_table_region(
region=list(range(len(mapping))), # reference row indices 0 to N-1
description="all electrodes",
)
from hdmf.backends.hdf5.h5_utils import H5DataIO
import numpy as np
import time
# Must transpose this data to be compatible with NWB
# Maxwell stores raw traces as (Nchans, Nsamples)
# NWB stores raw traces (Nsamples, Nchans)
raw = f['data_store/data0000/groups/routed/raw']
fname = '.tmp.bin'
tmpbin = np.memmap(datapath / fname,dtype=np.uint16,mode='w+',shape=(raw.shape[1],raw.shape[0]))
print("memmap formed.")
# Read/write in chunks
chunksize = 20000 * 10
Nchunks = raw.shape[1] // chunksize + 1 * (raw.shape[1] % chunksize != 0)
for chidx in range(Nchunks):
t0 = time.time()
idx0 = chidx * chunksize
idx1 = min((chidx+1) * chunksize, raw.shape[1])
X = np.array(raw[:,idx0:idx1],dtype=np.uint16)
tmpbin[idx0:idx1,:] = X.T
t1 = time.time()
print(f"Chunk {chidx+1} of {Nchunks} done in {t1-t0:.3f}")
# Write changes to disk
tmpbin.flush()
rawh5io = H5DataIO(data=tmpbin,compression='gzip',compression_opts=4,chunks=True)
raw_electrical_series = ElectricalSeries(
name="ElectricalSeries",
description="Raw acquisition traces",
data=rawh5io,
electrodes=all_table_region,
starting_time=0.0, # timestamp of the first sample in seconds relative to the session start time
rate=20000.0, # in Hz
)
nwbfile.add_acquisition(raw_electrical_series)
ecephys_module = nwbfile.create_processing_module(
name="ecephys", description="processed extracellular electrophysiology data"
)
from pynwb.ecephys import EventDetection
spks = f['data_store/data0000/spikes']
fr0 = f['data_store/data0000/groups/routed/frame_nos'][0]
thr = f['data_store/data0000/settings/spike_threshold'][0]
ch2idx = {ch: idx for idx,ch in enumerate(f['data_store/data0000/groups/routed/channels'])}
# 2 edge cases to consider. (i) non-routed channel detects a spike (ii) spikes are detected before frame 0
goodspks = [True if ch in ch2idx.keys() else False for ch in spks['channel']]
src_idxs = np.array([ch2idx[ch] for ch in spks['channel'][goodspks]])
src_idxs = np.column_stack([spks['frameno'][goodspks] - fr0,src_idxs])
src_idxs = src_idxs[np.argwhere(src_idxs[:,0] > 0)[0].squeeze():]
times = src_idxs[:,0]*50e-6
event_detection = EventDetection(
name="mxw_events",
detection_method=f"MaxWell built in RMS thresholding, {thr} * RMS",
source_electricalseries=raw_electrical_series,
source_idx=src_idxs,
times=times,
)
ecephys_module.add(event_detection)
with NWBHDF5IO("test_out.nwb", "w") as io:
io.write(nwbfile)
# With all of the above successful, we can clean up the tmp bin used to
# transpose the dataset
Path.unlink(datapath / fname)
Traceback
Validating div21-carbachol30um.nwb against cached namespace information using namespace 'core'.
- found the following errors:
EventDetection/source_idx (processing/ecephys/mxw_events/source_idx): incorrect shape - expected '[None]', got '(2469376, 2)'
Operating System
Linux
Python Executable
Python
Python Version
3.9
Package Versions
anyio @ file:///croot/anyio_1729121277521/work argon2-cffi @ file:///opt/conda/conda-bld/argon2-cffi_1645000214183/work argon2-cffi-bindings @ file:///tmp/build/80754af9/argon2-cffi-bindings_1644569679365/work arrow==1.3.0 asciitree==0.3.3 astroid @ file:///home/conda/feedstock_root/build_artifacts/astroid_1735074151157/work asttokens @ file:///opt/conda/conda-bld/asttokens_1646925590279/work async-lru @ file:///croot/async-lru_1699554519285/work attrs @ file:///croot/attrs_1729089401488/work autopep8 @ file:///home/conda/feedstock_root/build_artifacts/autopep8_1693061251004/work Babel @ file:///croot/babel_1671781930836/work backcall @ file:///home/ktietz/src/ci/backcall_1611930011877/work beautifulsoup4 @ file:///croot/beautifulsoup4-split_1718029820055/work bleach @ file:///opt/conda/conda-bld/bleach_1641577558959/work Bottleneck @ file:///croot/bottleneck_1731058641041/work Brotli @ file:///croot/brotli-split_1714483155106/work certifi @ file:///croot/certifi_1738623731865/work/certifi cffi @ file:///croot/cffi_1726856441404/work charset-normalizer @ file:///croot/charset-normalizer_1721748349566/work colorama @ file:///home/conda/feedstock_root/build_artifacts/colorama_1733218098505/work comm @ file:///croot/comm_1709322850197/work contourpy @ file:///croot/contourpy_1700583582875/work cuda-python==12.6.0 cycler @ file:///tmp/build/80754af9/cycler_1637851556182/work debugpy @ file:///croot/debugpy_1690905042057/work decorator @ file:///opt/conda/conda-bld/decorator_1643638310831/work defusedxml @ file:///tmp/build/80754af9/defusedxml_1615228127516/work dill @ file:///home/conda/feedstock_root/build_artifacts/dill_1733249551891/work distinctipy==1.3.4 docstring-to-markdown @ file:///home/conda/feedstock_root/build_artifacts/docstring-to-markdown_1733367229476/work ecos @ file:///croot/ecos_1677702164573/work elephant==1.1.1 entrypoints @ file:///home/conda/feedstock_root/build_artifacts/entrypoints_1643888246732/work exceptiongroup @ file:///croot/exceptiongroup_1706031385326/work executing @ file:///opt/conda/conda-bld/executing_1646925071911/work faiss-cpu==1.9.0 fasteners==0.19 fastjsonschema @ file:///opt/conda/conda-bld/python-fastjsonschema_1661371079312/work filelock==3.16.1 flake8 @ file:///home/conda/feedstock_root/build_artifacts/flake8_1733246708104/work fonttools @ file:///croot/fonttools_1713551344105/work fqdn==1.5.1 fsspec==2024.10.0 h11 @ file:///croot/h11_1706652277403/work h5py==3.13.0 hdmf==4.0.0 httpcore @ file:///croot/httpcore_1706728464539/work httpx @ file:///croot/httpx_1723474802858/work idna @ file:///croot/idna_1714398848350/work -e git+https://github.com/pwollstadt/IDTxl.git@8d3dacfabe92715e69032310dcb043d9a46d01b9#egg=idtxl importlib-metadata @ file:///croot/importlib_metadata-suite_1704813515092/work importlib_resources @ file:///croot/importlib_resources-suite_1720641103994/work ipykernel @ file:///croot/ipykernel_1737660677549/work ipympl @ file:///home/conda/feedstock_root/build_artifacts/ipympl_1713251546026/work ipython @ file:///croot/ipython_1694181358621/work ipython_genutils @ file:///home/conda/feedstock_root/build_artifacts/ipython_genutils_1716278396992/work ipywidgets @ file:///croot/ipywidgets_1709574692113/work isoduration==20.11.0 isort @ file:///home/conda/feedstock_root/build_artifacts/isort_1733236185843/work jedi @ file:///croot/jedi_1721058342488/work Jinja2 @ file:///croot/jinja2_1716993405101/work joblib==1.4.2 JPype1 @ file:///home/conda/feedstock_root/build_artifacts/jpype1_1704184898083/work json5 @ file:///croot/json5_1730786798687/work jsonpointer==3.0.0 jsonschema @ file:///croot/jsonschema_1728486696720/work jsonschema-specifications @ file:///croot/jsonschema-specifications_1699032386549/work jupyter @ file:///croot/jupyter_1707947101020/work jupyter-console @ file:///croot/jupyter_console_1679999630278/work jupyter-events @ file:///croot/jupyter_events_1718738097486/work jupyter-lsp @ file:///croot/jupyter-lsp-meta_1699978238815/work jupyter_client @ file:///croot/jupyter_client_1699455897726/work jupyter_core @ file:///croot/jupyter_core_1718818295206/work jupyter_server @ file:///croot/jupyter_server_1718827083372/work jupyter_server_terminals @ file:///croot/jupyter_server_terminals_1686870725608/work jupyterlab @ file:///home/conda/feedstock_root/build_artifacts/jupyterlab_1730308726474/work jupyterlab-lsp @ file:///home/conda/feedstock_root/build_artifacts/jupyter-lsp-meta_1696776403652/work/jupyterlab-lsp jupyterlab-pygments @ file:///tmp/build/80754af9/jupyterlab_pygments_1601490720602/work jupyterlab-vim @ file:///home/conda/feedstock_root/build_artifacts/jupyterlab_vim_1726493767827/work jupyterlab-widgets @ file:///croot/jupyterlab_widgets_1709322880313/work jupyterlab_server @ file:///croot/jupyterlab_server_1725865349919/work jupyternotify==0.1.15 kilosort==4.0.30 kiwisolver @ file:///croot/kiwisolver_1672387140495/work llvmlite==0.43.0 Mako @ file:///home/conda/feedstock_root/build_artifacts/mako_1715711344987/work MarkupSafe @ file:///croot/markupsafe_1704205993651/work matplotlib @ file:///tmp/build/80754af9/matplotlib-suite_1634667019719/work matplotlib-inline @ file:///opt/conda/conda-bld/matplotlib-inline_1662014470464/work mccabe @ file:///home/conda/feedstock_root/build_artifacts/mccabe_1733216466933/work mistune @ file:///opt/conda/conda-bld/mistune_1661496219659/work mpmath==1.3.0 nbclassic @ file:///home/conda/feedstock_root/build_artifacts/nbclassic_1716838762700/work nbclient @ file:///croot/nbclient_1698934205032/work nbconvert @ file:///home/conda/feedstock_root/build_artifacts/nbconvert-meta_1714477135335/work nbformat @ file:///croot/nbformat_1728049424075/work neo==0.13.4 nest-asyncio @ file:///croot/nest-asyncio_1708532673751/work networkx==3.2.1 notebook @ file:///home/conda/feedstock_root/build_artifacts/notebook_1680870634737/work notebook_shim @ file:///croot/notebook-shim_1699455894279/work numba==0.60.0 numcodecs==0.12.1 numexpr @ file:///croot/numexpr_1696515281613/work numpy==1.26.4 nvidia-cublas-cu12==12.4.5.8 nvidia-cuda-cupti-cu12==12.4.127 nvidia-cuda-nvrtc-cu12==12.4.127 nvidia-cuda-runtime-cu12==12.4.127 nvidia-cudnn-cu12==9.1.0.70 nvidia-cufft-cu12==11.2.1.3 nvidia-curand-cu12==10.3.5.147 nvidia-cusolver-cu12==11.6.1.9 nvidia-cusparse-cu12==12.3.1.170 nvidia-ml-py @ file:///home/conda/feedstock_root/build_artifacts/nvidia-ml-py_1734498123270/work nvidia-nccl-cu12==2.21.5 nvidia-nvjitlink-cu12==12.4.127 nvidia-nvtx-cu12==12.4.127 overrides @ file:///croot/overrides_1699371140756/work packaging @ file:///croot/packaging_1720101850331/work pandas==1.5.3 pandocfilters @ file:///opt/conda/conda-bld/pandocfilters_1643405455980/work parso @ file:///opt/conda/conda-bld/parso_1641458642106/work patsy @ file:///croot/patsy_1738159930729/work pexpect @ file:///tmp/build/80754af9/pexpect_1605563209008/work pickleshare @ file:///tmp/build/80754af9/pickleshare_1606932040724/work pillow @ file:///croot/pillow_1721059439630/work platformdirs @ file:///croot/platformdirs_1692205439124/work pluggy @ file:///home/conda/feedstock_root/build_artifacts/pluggy_1733222765875/work ply==3.11 probeinterface==0.2.24 prometheus-client @ file:///tmp/abs_d3zeliano1/croots/recipe/prometheus_client_1659455100375/work prompt-toolkit @ file:///croot/prompt-toolkit_1704404351921/work proplot @ file:///home/conda/feedstock_root/build_artifacts/proplot_1677584975715/work psutil @ file:///opt/conda/conda-bld/psutil_1656431268089/work ptyprocess @ file:///tmp/build/80754af9/ptyprocess_1609355006118/work/dist/ptyprocess-0.7.0-py2.py3-none-any.whl pure-eval @ file:///opt/conda/conda-bld/pure_eval_1646925070566/work py-cpuinfo @ file:///croot/py-cpuinfo_1698068081650/work pycodestyle @ file:///home/conda/feedstock_root/build_artifacts/pycodestyle_1733216196861/work pycparser @ file:///tmp/build/80754af9/pycparser_1636541352034/work pydocstyle @ file:///home/conda/feedstock_root/build_artifacts/pydocstyle_1733261631732/work pyflakes @ file:///home/conda/feedstock_root/build_artifacts/pyflakes_1733216066937/work Pygments @ file:///croot/pygments_1684279966437/work pyinform==0.2.0 pylint @ file:///home/conda/feedstock_root/build_artifacts/pylint_1735081221338/work pynvml @ file:///home/conda/feedstock_root/build_artifacts/pynvml_1734538935786/work pynwb==3.0.0 pyopencl @ file:///home/conda/feedstock_root/build_artifacts/pyopencl_1729563631734/work pyparsing @ file:///croot/pyparsing_1725041624202/work PyQt5==5.15.10 PyQt5-sip @ file:///croot/pyqt-split_1698769088074/work/pyqt_sip PyQt6==6.7.1 PyQt6-Qt6==6.7.3 PyQt6_sip==13.8.0 pyqtgraph==0.13.7 PySocks @ file:///tmp/build/80754af9/pysocks_1605305812635/work python-dateutil @ file:///croot/python-dateutil_1716495738603/work python-json-logger @ file:///croot/python-json-logger_1683823803357/work python-lsp-jsonrpc @ file:///home/conda/feedstock_root/build_artifacts/python-lsp-jsonrpc_1736013464386/work python-lsp-server @ file:///home/conda/feedstock_root/build_artifacts/python-lsp-server-meta_1736028135583/work pytoolconfig @ file:///home/conda/feedstock_root/build_artifacts/pytoolconfig_1733360721313/work pytools @ file:///home/conda/feedstock_root/build_artifacts/pytools_1725633725292/work pytz @ file:///croot/pytz_1713974312559/work PyYAML @ file:///croot/pyyaml_1728657952215/work pyzmq @ file:///croot/pyzmq_1705605076900/work qtconsole @ file:///croot/qtconsole_1726836502164/work QtPy @ file:///croot/qtpy_1700144840038/work quantities==0.16.1 referencing @ file:///croot/referencing_1699012038513/work requests @ file:///croot/requests_1721410876868/work rfc3339-validator @ file:///croot/rfc3339-validator_1683077044675/work rfc3986-validator @ file:///croot/rfc3986-validator_1683058983515/work -e git+ssh://[email protected]/Kasun-Pathirage/rol_project.git@7f157b08d0c58e08d1e07ec1f666210140985fd4#egg=rolpy rope @ file:///home/conda/feedstock_root/build_artifacts/rope_1733342323938/work rpds-py @ file:///croot/rpds-py_1698945930462/work ruamel.yaml==0.18.10 ruamel.yaml.clib==0.2.12 scikit-learn==1.5.2 scipy @ file:///croot/scipy_1701295040508/work/dist/scipy-1.11.4-cp39-cp39-linux_x86_64.whl#sha256=7d3c50abd59d682a73a495b68c974347054667b7f0eee3dcdf2605778b9e5c6f Send2Trash @ file:///croot/send2trash_1699371139552/work sip @ file:///croot/sip_1698675935381/work siphash24 @ file:///home/conda/feedstock_root/build_artifacts/siphash24_1729087984596/work/dist/siphash24-1.7-cp39-cp39-linux_x86_64.whl#sha256=7f64e041d78391953c05644e7ced6e022e7a5fa999a4187428dea9a6fd63551d six @ file:///tmp/build/80754af9/six_1644875935023/work sniffio @ file:///croot/sniffio_1705431295498/work snowballstemmer @ file:///home/conda/feedstock_root/build_artifacts/snowballstemmer_1637143057757/work soupsieve @ file:///croot/soupsieve_1696347547217/work spikeinterface==0.101.2 spycon @ git+https://github.com/christiando/spycon@dc9afb6c9b5f030d51a6fd70cff56f3d8db59983 stack-data @ file:///opt/conda/conda-bld/stack_data_1646927590127/work statsmodels @ file:///croot/statsmodels_1689937266057/work sympy==1.13.1 tables @ file:///croot/pytables_1705614842881/work terminado @ file:///croot/terminado_1671751832461/work threadpoolctl==3.5.0 tinycss2 @ file:///croot/tinycss2_1668168815555/work tomli @ file:///opt/conda/conda-bld/tomli_1657175507142/work tomlkit @ file:///home/conda/feedstock_root/build_artifacts/tomlkit_1733230743009/work torch==2.5.0 tornado @ file:///croot/tornado_1718740109488/work tqdm==4.66.5 traitlets @ file:///home/conda/feedstock_root/build_artifacts/traitlets_1675110562325/work triton==3.1.0 types-python-dateutil==2.9.0.20241206 typing_extensions @ file:///croot/typing_extensions_1715268824938/work tzdata==2024.2 ujson @ file:///home/conda/feedstock_root/build_artifacts/ujson_1724954423010/work unicodedata2 @ file:///croot/unicodedata2_1713212950228/work uri-template==1.3.0 urllib3 @ file:///croot/urllib3_1727769808118/work wcwidth @ file:///Users/ktietz/demo/mc3/conda-bld/wcwidth_1629357192024/work webcolors==24.11.1 webencodings==0.5.1 websocket-client @ file:///croot/websocket-client_1715878298792/work whatthepatch @ file:///home/conda/feedstock_root/build_artifacts/whatthepatch_1736040839611/work widgetsnbextension @ file:///croot/widgetsnbextension_1709322880396/work xarray @ file:///croot/xarray_1689041471350/work yapf @ file:///home/conda/feedstock_root/build_artifacts/yapf_1736192735449/work zarr==2.17.2 zipp @ file:///croot/zipp_1729012354496/work
Code of Conduct
- [x] I agree to follow this project's Code of Conduct
- [x] Have you checked the Contributing document?
- [x] Have you ensured this bug was not already reported?
Hi @naterenegar, thanks for the detailed issue!
I'm attempting to convert my data to the NWB format, and have some voltage traces recorded from a multi-electrode array. Spikes were detected in these traces. If I want to save the traces and the spike events, the way to go seems like I should use the EventDetection data interface. If my raw data is two-dimensional, e.g., 100,000 samples and 1000 channels, would the source_idx also be two dimensional?
You are correct that EventDetection is a good way to save information about when spikes were detected from your raw ElectricalSeries data. Unfortunately the current definition of EventDetection only supports specifying the index along the time axis of an event, and not the specific channel on which the event was detected. This is something we should update in the schema and have created an issue and related PR for: https://github.com/NeurodataWithoutBorders/nwb-schema/issues/619.
If you need a temporary (though not ideal) workaround, there are a couple of options:
- You can make two
EventDetectionobjects where one uses thesource_indexto indicate the sample index along the time dimension and another uses thesource_indexto indicate the channel index. - You can use the
SpikeEventSeriesdata type which can have the dimensions [num_events, num_channels, num_samples] as a sparse array, where the channel on which the event was detected has a single sample value of 1 (usuallySpikeEventSeriesare used to store the spike traces for each event detected across channels).
Support for EventDetection with a two-dimensional source_idx shape of (N, 2) was added in #2091. These updates have now been merged into the dev branch and will be included in the upcoming pynwb release.
I am closing this issue for now since these updates should support the original use case described here, but please reopen if you have any further questions!