From 720eb00e60c6df928f4e6cbe938d4db625feab58 Mon Sep 17 00:00:00 2001 From: Chris Broz Date: Wed, 10 Nov 2021 14:26:24 -0600 Subject: [PATCH 01/77] Mult rootdirs. If sess_dir, check fullpath. Give OpenEphys fullpath. --- element_array_ephys/ephys.py | 17 +++++++++++------ 1 file changed, 11 insertions(+), 6 deletions(-) diff --git a/element_array_ephys/ephys.py b/element_array_ephys/ephys.py index 3eec7842..1cd2ee00 100644 --- a/element_array_ephys/ephys.py +++ b/element_array_ephys/ephys.py @@ -141,20 +141,22 @@ class EphysFile(dj.Part): def make(self, key): sess_dir = pathlib.Path(get_session_directory(key)) + sess_dir_full = find_full_path(get_ephys_root_data_dir(), sess_dir) inserted_probe_serial_number = (ProbeInsertion * probe.Probe & key).fetch1('probe') # search session dir and determine acquisition software for ephys_pattern, ephys_acq_type in zip(['*.ap.meta', '*.oebin'], ['SpikeGLX', 'Open Ephys']): - ephys_meta_filepaths = [fp for fp in sess_dir.rglob(ephys_pattern)] + ephys_meta_filepaths = [fp for fp in sess_dir_full.rglob(ephys_pattern)] if ephys_meta_filepaths: acq_software = ephys_acq_type break else: raise FileNotFoundError( f'Ephys recording data not found!' - f' Neither SpikeGLX nor Open Ephys recording files found') + f' Neither SpikeGLX nor Open Ephys recording files found' + f' in {sess_dir}') if acq_software == 'SpikeGLX': for meta_filepath in ephys_meta_filepaths: @@ -192,7 +194,7 @@ def make(self, key): **key, 'file_path': meta_filepath.relative_to(root_dir).as_posix()}) elif acq_software == 'Open Ephys': - dataset = openephys.OpenEphys(sess_dir) + dataset = openephys.OpenEphys(sess_dir_full) for serial_number, probe_data in dataset.probes.items(): if str(serial_number) == inserted_probe_serial_number: break @@ -291,7 +293,8 @@ def make(self, key): electrode_keys.append(probe_electrodes[(shank, shank_col, shank_row)]) elif acq_software == 'Open Ephys': sess_dir = pathlib.Path(get_session_directory(key)) - loaded_oe = openephys.OpenEphys(sess_dir) + sess_dir_full = find_full_path(get_ephys_root_data_dir(), sess_dir) + loaded_oe = openephys.OpenEphys(sess_dir_full) oe_probe = loaded_oe.probes[probe_sn] lfp_channel_ind = np.arange( @@ -614,7 +617,8 @@ def yield_unit_waveforms(): neuropixels_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) elif acq_software == 'Open Ephys': sess_dir = pathlib.Path(get_session_directory(key)) - openephys_dataset = openephys.OpenEphys(sess_dir) + sess_dir_full = find_full_path(get_ephys_root_data_dir(), sess_dir) + openephys_dataset = openephys.OpenEphys(sess_dir_full) neuropixels_recording = openephys_dataset.probes[probe_serial_number] def yield_unit_waveforms(): @@ -697,7 +701,8 @@ def get_neuropixels_channel2electrode_map(ephys_recording_key, acq_software): spikeglx_meta.shankmap['data'])} elif acq_software == 'Open Ephys': sess_dir = pathlib.Path(get_session_directory(ephys_recording_key)) - openephys_dataset = openephys.OpenEphys(sess_dir) + sess_dir_full = find_full_path(get_ephys_root_data_dir(), sess_dir) + openephys_dataset = openephys.OpenEphys(sess_dir_full) probe_serial_number = (ProbeInsertion & ephys_recording_key).fetch1('probe') probe_dataset = openephys_dataset.probes[probe_serial_number] From f3dd8d80c0a0b9a4e16734a35cd0e2f4520a0142 Mon Sep 17 00:00:00 2001 From: bendichter Date: Wed, 15 Dec 2021 14:45:47 -0500 Subject: [PATCH 02/77] add draft convert to nwb --- element_array_ephys/export/__init__.py | 0 element_array_ephys/export/nwb.py | 384 +++++++++++++++++++++++++ 2 files changed, 384 insertions(+) create mode 100644 element_array_ephys/export/__init__.py create mode 100644 element_array_ephys/export/nwb.py diff --git a/element_array_ephys/export/__init__.py b/element_array_ephys/export/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py new file mode 100644 index 00000000..064a0501 --- /dev/null +++ b/element_array_ephys/export/nwb.py @@ -0,0 +1,384 @@ +from hdmf.backends.hdf5 import H5DataIO +from nwb_conversion_tools.utils.genericdatachunkiterator import GenericDataChunkIterator +from nwb_conversion_tools.utils.recordingextractordatachunkiterator import ( + RecordingExtractorDataChunkIterator, +) +from pynwb import NWBHDF5IO +from pynwb.ecephys import ElectricalSeries, LFP +from pynwb.misc import Units +from spikeextractors import ( + SpikeGLXRecordingExtractor, + SubRecordingExtractor, + OpenEphysNPIXRecordingExtractor, +) +from tqdm import tqdm + +from element_session import session + + +# from element_data_loader.utils import find_full_path + + +class LFPDataChunkIterator(GenericDataChunkIterator): + """DataChunkIterator for LFP data that pulls data one channel at a time.""" + + def __init__( + self, + lfp_electrodes_query, + buffer_length: int = None, + chunk_length: int = 10000, + ): + self.lfp_electrodes_query = lfp_electrodes_query + self.electrodes = self.lfp_electrodes_query.fetch("electrode") + + first_record = ( + self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) + ).fetch(as_dict=True)[0] + + self.n_channels = len(lfp_electrodes_query) + self.n_tt = len(first_record["lfp"]) + self._dtype = first_record["lfp"].dtype + + chunk_shape = (chunk_length, 1) + if buffer_length is not None: + buffer_shape = (buffer_length, 1) + else: + buffer_shape = (len(first_record["lfp"]), 1) + super().__init__(buffer_shape=buffer_shape, chunk_shape=chunk_shape) + + def _get_data(self, selection): + + electrode = self.electrodes[selection[1]][0] + return (self.lfp_electrodes_query & dict(electrode=electrode)).fetch1("lfp")[ + selection[0], np.newaxis + ] + + def _get_dtype(self): + return self._dtype + + def _get_maxshape(self): + return self.n_tt, self.n_channels + + +def check_module(nwbfile, name: str, description: str = None): + """ + Check if processing module exists. If not, create it. Then return module. + Parameters + ---------- + nwbfile: pynwb.NWBFile + name: str + description: str | None (optional) + Returns + ------- + pynwb.module + """ + assert isinstance(nwbfile, NWBFile), "'nwbfile' should be of type pynwb.NWBFile" + if name in nwbfile.processing: + return nwbfile.processing[name] + else: + if description is None: + description = name + return nwbfile.create_processing_module(name, description) + + +def get_ephys_root_data_dir(): + root_data_dirs = dj.config.get("custom", {}).get("ephys_root_data_dir", None) + return root_data_dirs + + +def add_electrodes_to_nwb(session_key: dict, nwbfile: NWBFile): + + electrodes_query = probe.ProbeType.Electrode * probe.ElectrodeConfig.Electrode + + session_electrodes = ( + probe.ElectrodeConfig.Electrode & (ephys.EphysRecording & session_key) + ).fetch() + + for additional_attribute in ["shank_col", "shank_row", "shank"]: + nwbfile.add_electrode_column( + name=electrodes_query.heading.attributes[additional_attribute].name, + description=electrodes_query.heading.attributes[ + additional_attribute + ].comment, + ) + + nwbfile.add_electrode_column( + name="id_in_probe", description="electrode id within the probe", + ) + + session_probes = (ephys.ProbeInsertion * probe.Probe & session_key).fetch( + as_dict=True + ) + for this_probe in session_probes: + + insertion_record = (ephys.InsertionLocation & this_probe).fetch() + if False: # insertion_record: + insert_location = insertion_record.fetch(*non_primary_keys, as_dict=True) + insert_location = json.dumps(insert_location) + else: + insert_location = "unknown" + + device = nwbfile.create_device( + name=this_probe["probe"], + description=this_probe.get("probe_comment", None), + manufacturer=this_probe.get("probe_type", None), + ) + shank_ids = np.unique((probe.ProbeType.Electrode & this_probe).fetch("shank")) + for shank_id in shank_ids: + electrode_group = nwbfile.create_electrode_group( + name=f"probe{this_probe['probe']}_shank{shank_id}", + description=f"probe{this_probe['probe']}_shank{shank_id}", + location=insert_location, + device=device, + ) + + electrodes_query = ( + probe.ProbeType.Electrode & this_probe & dict(shank=shank_id) + ).fetch(as_dict=True) + for electrode in electrodes_query: + nwbfile.add_electrode( + group=electrode_group, + filtering="unknown", + imp=-1.0, + x=np.nan, # to do: populate these values once the CCF element is ready + y=np.nan, + z=np.nan, + rel_x=electrode["x_coord"], + rel_y=electrode["y_coord"], + rel_z=np.nan, + shank_col=electrode["shank_col"], + shank_row=electrode["shank_row"], + location="unknown", + id_in_probe=electrode["electrode"], + shank=electrode["shank"], + ) + + +def create_units_table( + session_key, + nwbfile, + units_query, + paramset_record, + name="units", + desc="data on spiking units", +): + + # electrode id mapping + electrode_id_mapping = { + nwbfile.electrodes.id.data[x]: x for x in range(len(nwbfile.electrodes.id)) + } + + units_table = Units(name=name, description=desc) + for additional_attribute in ["cluster_quality_label", "spike_depths"]: + units_table.add_column( + name=units_query.heading.attributes[additional_attribute].name, + description=units_query.heading.attributes[additional_attribute].comment, + index=True if additional_attribute == "spike_depths" else False, + ) + + clustering_query = ( + ephys.EphysRecording * ephys.ClusteringTask & session_key & paramset_record + ) + + for unit in tqdm( + (clustering_query @ ephys.CuratedClustering.Unit).fetch(as_dict=True), + desc=f"creating units table for paramset {paramset_record['paramset_idx']}", + ): + + probe_id, shank_num = ( + probe.ProbeType.Electrode + * ephys.ProbeInsertion + * ephys.CuratedClustering.Unit + & {"unit": unit["unit"], "insertion_number": unit["insertion_number"]} + ).fetch1("probe", "shank") + + waveform_mean = ( + ephys.WaveformSet.PeakWaveform() & clustering_query & unit + ).fetch1("peak_electrode_waveform") + + units_table.add_row( + id=unit["unit"], + electrodes=[electrode_id_mapping[unit["electrode"]]], + electrode_group=nwbfile.electrode_groups[ + f"probe{probe_id}_shank{shank_num}" + ], + cluster_quality_label=unit["cluster_quality_label"], + spike_times=unit["spike_times"], + spike_depths=unit["spike_depths"], + waveform_mean=waveform_mean, + ) + + return units_table + + +def add_ephys_units_to_nwb( + session_key: dict, nwbfile: NWBFile, primary_clustering_paramset_idx: int = 0 +): + + if not ephys.ClusteringTask & session_key: + return + + if nwbfile.electrodes is None: + add_electrodes_to_nwb(session_key, nwbfile) + + # add additional columns to the units table + units_query = ephys.CuratedClustering.Unit() & session_key + + for paramset_record in ( + ephys.ClusteringParamSet & ephys.CuratedClustering & session_key + ).fetch("paramset_idx", "clustering_method", "paramset_desc", as_dict=True): + if paramset_record["paramset_idx"] == primary_clustering_paramset_idx: + units_table = create_units_table( + session_key, + nwbfile, + units_query, + paramset_record, + desc=paramset_record["paramset_desc"], + ) + nwbfile.units = units_table + else: + name = f"units_{paramset_record['clustering_method']}" + units_table = create_units_table( + session_key, + units_query, + paramset_record, + name=name, + desc=paramset_record["paramset_desc"], + ) + ecephys_module = check_module(nwbfile, "ecephys") + ecephys_module.add(units_table) + + +def add_ephys_recording_to_nwb( + session_key: dict, nwbfile: NWBFile, end_frame: int = None +): + + if not ephys.EphysRecording & session_key: + return + + if nwbfile.electrodes is None: + add_electrodes_to_nwb(session_key, nwbfile) + + mapping = { + ( + nwbfile.electrodes["group"][idx].device.name, + nwbfile.electrodes["id_in_probe"][idx], + ): idx + for idx in range(len(nwbfile.electrodes)) + } + + for ephys_recording_record in (ephys.EphysRecording & session_key).fetch( + as_dict=True + ): + probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") + + # relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") + # file_path = find_full_path(get_ephys_root_data_dir(), relative_path) + file_path = "../inbox/subject5/session1/probe_1/npx_g0_t0.imec.ap.bin" + + if ephys_recording_record["acq_software"] == "SpikeGLX": + extractor = SpikeGLXRecordingExtractor(file_path=file_path) + elif ephys_recording_record["acq_software"] == "OpenEphys": + extractor = OpenEphysNPIXRecordingExtractor(file_path=file_path) + channel_conversions = extractor.get_channel_gains() + # to do: channel conversions for OpenEphys + + if end_frame is not None: + extractor = SubRecordingExtractor(extractor, end_frame=end_frame) + + recording_channels_by_id = ( + probe.ElectrodeConfig.Electrode() + & ephys.EphysRecording() + & session_key + & dict(insertion_number=1) + ).fetch("electrode") + + nwbfile.add_acquisition( + ElectricalSeries( # to do: add conversion + name=f"ElectricalSeries{ephys_recording_record['insertion_number']}", + description=str(ephys_recording_record), + data=H5DataIO( + RecordingExtractorDataChunkIterator(extractor), compression=True + ), + rate=ephys_recording_record["sampling_rate"], + starting_time=( + ephys_recording_record["recording_datetime"] + - ephys_recording_record["session_datetime"] + ).total_seconds(), + electrodes=nwbfile.create_electrode_table_region( + region=[mapping[(probe_id, x)] for x in recording_channels_by_id], + name="electrodes", + description="recorded electrodes", + ), + conversion=1e-6, + channel_conversion=channel_conversions, + ) + ) + + +def add_ephys_lfp_to_nwb(session_key: dict, nwbfile: NWBFile): + + if not ephys.LFP & session_key: + return + + if nwbfile.electrodes is None: + add_electrodes_to_nwb(session_key, nwbfile) + + ecephys_module = check_module( + nwbfile, name="ecephys", description="preprocessed ephys data" + ) + + nwb_lfp = LFP(name="LFP") + ecephys_module.add(nwb_lfp) + + mapping = { + ( + nwbfile.electrodes["group"][idx].device.name, + nwbfile.electrodes["id_in_probe"][idx], + ): idx + for idx in range(len(nwbfile.electrodes)) + } + + for lfp_record in (ephys.LFP & session_key).fetch(as_dict=True): + probe_id = (ephys.ProbeInsertion & lfp_record).fetch1("probe") + + lfp_electrodes_query = ephys.LFP.Electrode & lfp_record + lfp_data = LFPDataChunkIterator(lfp_electrodes_query) + + nwb_lfp.create_electrical_series( + name=f"ElectricalSeries{lfp_record['insertion_number']}", + description=f"LFP from probe {probe_id}", + data=H5DataIO(lfp_data, compression=True), + timestamps=lfp_record["lfp_time_stamps"], + electrodes=nwbfile.create_electrode_table_region( + name="electrodes", + description="electrodes used for LFP", + region=[ + mapping[(probe_id, x)] + for x in lfp_electrodes_query.fetch("electrode") + ], + ), + ) + +def session_to_nwb(session_key, subject_id=None, raw=True, spikes=True, lfp=True, end_frame=None): + + nwbfile = session.export.nwb.session_to_nwb(session_key, subject_id=subject_id) + + if raw: + add_ephys_recording_to_nwb(session_key, nwbfile, end_frame=end_frame) + + if spikes: + add_ephys_units_to_nwb(session_key, nwbfile) + + if lfp: + add_ephys_lfp_to_nwb(session_key, nwbfile) + + return nwbfile + +def write_nwb(nwbfile, fname, check_read=True): + with NWBHDF5IO(fname, 'w') as io: + io.write(nwbfile) + + if check_read: + with NWBHDF5IO(fname, 'r') as io: + io.read() From 8d75a61ac29b0d51381104472a8570db6c99bcb5 Mon Sep 17 00:00:00 2001 From: bendichter Date: Wed, 5 Jan 2022 16:09:06 -0500 Subject: [PATCH 03/77] * add lfp from source * docstrings * json dump insertion into location * ignore channel_conversion if all are 1 * black formatting --- element_array_ephys/export/nwb.py | 369 +++++++++++++++++++++--------- 1 file changed, 265 insertions(+), 104 deletions(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index 064a0501..dd85e8c8 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -1,11 +1,13 @@ +import decimal +import numpy as np from hdmf.backends.hdf5 import H5DataIO from nwb_conversion_tools.utils.genericdatachunkiterator import GenericDataChunkIterator from nwb_conversion_tools.utils.recordingextractordatachunkiterator import ( RecordingExtractorDataChunkIterator, ) -from pynwb import NWBHDF5IO -from pynwb.ecephys import ElectricalSeries, LFP -from pynwb.misc import Units +from nwb_conversion_tools.utils.spike_interface import check_module + +import pynwb from spikeextractors import ( SpikeGLXRecordingExtractor, SubRecordingExtractor, @@ -19,32 +21,39 @@ # from element_data_loader.utils import find_full_path +class DecimalEncoder(json.JSONEncoder): + def default(self, o): + if isinstance(o, decimal.Decimal): + return str(o) + return super(DecimalEncoder, self).default(o) + + class LFPDataChunkIterator(GenericDataChunkIterator): - """DataChunkIterator for LFP data that pulls data one channel at a time.""" + """ + DataChunkIterator for LFP data that pulls data one channel at a time. Used when + reading LFP data from the database (as opposed to directly from source files) + """ - def __init__( - self, - lfp_electrodes_query, - buffer_length: int = None, - chunk_length: int = 10000, - ): + def __init__(self, lfp_electrodes_query, chunk_length: int = 10000): + """ + Parameters + ---------- + lfp_electrodes_query: element_array_ephys.ephys_no_curation.LFP + chunk_length: int, optional + Chunks are blocks of disk space where data are stored contiguously and compressed + """ self.lfp_electrodes_query = lfp_electrodes_query self.electrodes = self.lfp_electrodes_query.fetch("electrode") first_record = ( - self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) + self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) ).fetch(as_dict=True)[0] self.n_channels = len(lfp_electrodes_query) self.n_tt = len(first_record["lfp"]) self._dtype = first_record["lfp"].dtype - chunk_shape = (chunk_length, 1) - if buffer_length is not None: - buffer_shape = (buffer_length, 1) - else: - buffer_shape = (len(first_record["lfp"]), 1) - super().__init__(buffer_shape=buffer_shape, chunk_shape=chunk_shape) + super().__init__(buffer_shape=(self.n_tt, 1), chunk_shape=(chunk_length, 1)) def _get_data(self, selection): @@ -60,38 +69,21 @@ def _get_maxshape(self): return self.n_tt, self.n_channels -def check_module(nwbfile, name: str, description: str = None): +def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): """ - Check if processing module exists. If not, create it. Then return module. + Add electrodes table to NWBFile. This is needed for any ElectricalSeries, including + raw source data and LFP. + Parameters ---------- + session_key: dict nwbfile: pynwb.NWBFile - name: str - description: str | None (optional) - Returns - ------- - pynwb.module """ - assert isinstance(nwbfile, NWBFile), "'nwbfile' should be of type pynwb.NWBFile" - if name in nwbfile.processing: - return nwbfile.processing[name] - else: - if description is None: - description = name - return nwbfile.create_processing_module(name, description) - - -def get_ephys_root_data_dir(): - root_data_dirs = dj.config.get("custom", {}).get("ephys_root_data_dir", None) - return root_data_dirs - - -def add_electrodes_to_nwb(session_key: dict, nwbfile: NWBFile): electrodes_query = probe.ProbeType.Electrode * probe.ElectrodeConfig.Electrode session_electrodes = ( - probe.ElectrodeConfig.Electrode & (ephys.EphysRecording & session_key) + probe.ElectrodeConfig.Electrode & (ephys.EphysRecording & session_key) ).fetch() for additional_attribute in ["shank_col", "shank_row", "shank"]: @@ -106,15 +98,13 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: NWBFile): name="id_in_probe", description="electrode id within the probe", ) - session_probes = (ephys.ProbeInsertion * probe.Probe & session_key).fetch( - as_dict=True - ) - for this_probe in session_probes: - - insertion_record = (ephys.InsertionLocation & this_probe).fetch() - if False: # insertion_record: - insert_location = insertion_record.fetch(*non_primary_keys, as_dict=True) - insert_location = json.dumps(insert_location) + for this_probe in (ephys.ProbeInsertion * probe.Probe & session_key).fetch( + as_dict=True + ): + insertion_record = (ephys.InsertionLocation & this_probe).fetch1() + if insertion_record: + [insertion_record.pop(k) for k in ephys.InsertionLocation.heading.primary_key] + insert_location = json.dumps(insertion_record, cls=DecimalEncoder) else: insert_location = "unknown" @@ -133,7 +123,7 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: NWBFile): ) electrodes_query = ( - probe.ProbeType.Electrode & this_probe & dict(shank=shank_id) + probe.ProbeType.Electrode & this_probe & dict(shank=shank_id) ).fetch(as_dict=True) for electrode in electrodes_query: nwbfile.add_electrode( @@ -155,20 +145,31 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: NWBFile): def create_units_table( - session_key, - nwbfile, - units_query, - paramset_record, - name="units", - desc="data on spiking units", + session_key: dict, + nwbfile: pynwb.NWBFile, + units_query: ephys.CuratedClustering.Unit, + paramset_record, + name="units", + desc="data on spiking units", ): + """ + + Parameters + ---------- + session_key: dict + nwbfile: pynwb.NWBFile + units_query: ephys.CuratedClustering.Unit + paramset_record: int + name: str, optional + default="units" + desc: str, optional + default="data on spiking units" + """ # electrode id mapping - electrode_id_mapping = { - nwbfile.electrodes.id.data[x]: x for x in range(len(nwbfile.electrodes.id)) - } + mapping = get_electrodes_mapping(nwbfile.electrodes) - units_table = Units(name=name, description=desc) + units_table = pynwb.misc.Units(name=name, description=desc) for additional_attribute in ["cluster_quality_label", "spike_depths"]: units_table.add_column( name=units_query.heading.attributes[additional_attribute].name, @@ -177,28 +178,28 @@ def create_units_table( ) clustering_query = ( - ephys.EphysRecording * ephys.ClusteringTask & session_key & paramset_record + ephys.EphysRecording * ephys.ClusteringTask & session_key & paramset_record ) for unit in tqdm( - (clustering_query @ ephys.CuratedClustering.Unit).fetch(as_dict=True), - desc=f"creating units table for paramset {paramset_record['paramset_idx']}", + (clustering_query @ ephys.CuratedClustering.Unit).fetch(as_dict=True), + desc=f"creating units table for paramset {paramset_record['paramset_idx']}", ): probe_id, shank_num = ( - probe.ProbeType.Electrode - * ephys.ProbeInsertion - * ephys.CuratedClustering.Unit - & {"unit": unit["unit"], "insertion_number": unit["insertion_number"]} + probe.ProbeType.Electrode + * ephys.ProbeInsertion + * ephys.CuratedClustering.Unit + & {"unit": unit["unit"], "insertion_number": unit["insertion_number"]} ).fetch1("probe", "shank") waveform_mean = ( - ephys.WaveformSet.PeakWaveform() & clustering_query & unit + ephys.WaveformSet.PeakWaveform() & clustering_query & unit ).fetch1("peak_electrode_waveform") units_table.add_row( id=unit["unit"], - electrodes=[electrode_id_mapping[unit["electrode"]]], + electrodes=[mapping[(probe_id, unit["electrode"])]], electrode_group=nwbfile.electrode_groups[ f"probe{probe_id}_shank{shank_num}" ], @@ -212,8 +213,19 @@ def create_units_table( def add_ephys_units_to_nwb( - session_key: dict, nwbfile: NWBFile, primary_clustering_paramset_idx: int = 0 + session_key: dict, nwbfile: pynwb.NWBFile, primary_clustering_paramset_idx: int = 0 ): + """ + Add spiking data to NWBFile. + + Parameters + ---------- + session_key: dict + nwbfile: pynwb.NWBFile + primary_clustering_paramset_idx: int, optional + + + """ if not ephys.ClusteringTask & session_key: return @@ -225,7 +237,7 @@ def add_ephys_units_to_nwb( units_query = ephys.CuratedClustering.Unit() & session_key for paramset_record in ( - ephys.ClusteringParamSet & ephys.CuratedClustering & session_key + ephys.ClusteringParamSet & ephys.CuratedClustering & session_key ).fetch("paramset_idx", "clustering_method", "paramset_desc", as_dict=True): if paramset_record["paramset_idx"] == primary_clustering_paramset_idx: units_table = create_units_table( @@ -249,9 +261,39 @@ def add_ephys_units_to_nwb( ecephys_module.add(units_table) +def get_electrodes_mapping(electrodes): + """ + Create a mapping from the group (shank) and electrode id within that group to the row number of the electrodes + table. This is used in the construction of the DynamicTableRegion that indicates what rows of the electrodes + table correspond to the data in an ElectricalSeries. + + Parameters + ---------- + electrodes: hdmf.common.table.DynamicTable + + Returns + ------- + dict + + """ + return { + (electrodes["group"][idx].device.name, electrodes["id_in_probe"][idx],): idx + for idx in range(len(electrodes)) + } + + def add_ephys_recording_to_nwb( - session_key: dict, nwbfile: NWBFile, end_frame: int = None + session_key: dict, nwbfile: pynwb.NWBFile, end_frame: int = None ): + """ + + Parameters + ---------- + session_key: dict + nwbfile: NWBFile + end_frame: int, optional + Used for small test conversions + """ if not ephys.EphysRecording & session_key: return @@ -259,16 +301,10 @@ def add_ephys_recording_to_nwb( if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) - mapping = { - ( - nwbfile.electrodes["group"][idx].device.name, - nwbfile.electrodes["id_in_probe"][idx], - ): idx - for idx in range(len(nwbfile.electrodes)) - } + mapping = get_electrodes_mapping(nwbfile.electrodes) for ephys_recording_record in (ephys.EphysRecording & session_key).fetch( - as_dict=True + as_dict=True ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") @@ -280,21 +316,25 @@ def add_ephys_recording_to_nwb( extractor = SpikeGLXRecordingExtractor(file_path=file_path) elif ephys_recording_record["acq_software"] == "OpenEphys": extractor = OpenEphysNPIXRecordingExtractor(file_path=file_path) - channel_conversions = extractor.get_channel_gains() - # to do: channel conversions for OpenEphys + else: + raise ValueError( + f"unsupported acq_software type: {ephys_recording_record['acq_software']}" + ) + + if all(extractor.get_channel_gains() == 1): + channel_conversion = None + else: + channel_conversion = extractor.get_channel_gains() if end_frame is not None: extractor = SubRecordingExtractor(extractor, end_frame=end_frame) recording_channels_by_id = ( - probe.ElectrodeConfig.Electrode() - & ephys.EphysRecording() - & session_key - & dict(insertion_number=1) + probe.ElectrodeConfig.Electrode() & ephys_recording_record ).fetch("electrode") nwbfile.add_acquisition( - ElectricalSeries( # to do: add conversion + pynwb.ecephys.ElectricalSeries( name=f"ElectricalSeries{ephys_recording_record['insertion_number']}", description=str(ephys_recording_record), data=H5DataIO( @@ -302,8 +342,8 @@ def add_ephys_recording_to_nwb( ), rate=ephys_recording_record["sampling_rate"], starting_time=( - ephys_recording_record["recording_datetime"] - - ephys_recording_record["session_datetime"] + ephys_recording_record["recording_datetime"] + - ephys_recording_record["session_datetime"] ).total_seconds(), electrodes=nwbfile.create_electrode_table_region( region=[mapping[(probe_id, x)] for x in recording_channels_by_id], @@ -311,12 +351,20 @@ def add_ephys_recording_to_nwb( description="recorded electrodes", ), conversion=1e-6, - channel_conversion=channel_conversions, + channel_conversion=channel_conversion, ) ) -def add_ephys_lfp_to_nwb(session_key: dict, nwbfile: NWBFile): +def add_ephys_lfp_from_dj_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): + """ + Read LFP data from the data in element-aray-ephys + + Parameters + ---------- + session_key: dict + nwbfile: NWBFile + """ if not ephys.LFP & session_key: return @@ -328,16 +376,10 @@ def add_ephys_lfp_to_nwb(session_key: dict, nwbfile: NWBFile): nwbfile, name="ecephys", description="preprocessed ephys data" ) - nwb_lfp = LFP(name="LFP") + nwb_lfp = pynwb.ecephys.LFP(name="LFP") ecephys_module.add(nwb_lfp) - mapping = { - ( - nwbfile.electrodes["group"][idx].device.name, - nwbfile.electrodes["id_in_probe"][idx], - ): idx - for idx in range(len(nwbfile.electrodes)) - } + mapping = get_electrodes_mapping(nwbfile.electrodes) for lfp_record in (ephys.LFP & session_key).fetch(as_dict=True): probe_id = (ephys.ProbeInsertion & lfp_record).fetch1("probe") @@ -360,7 +402,111 @@ def add_ephys_lfp_to_nwb(session_key: dict, nwbfile: NWBFile): ), ) -def session_to_nwb(session_key, subject_id=None, raw=True, spikes=True, lfp=True, end_frame=None): + +def add_ephys_lfp_from_source_to_nwb( + session_key: dict, nwbfile: pynwb.NWBFile, end_frame=None +): + """ + Read the LFP data directly from the source file. Currently only works for SpikeGLX data. + + Parameters + ---------- + session_key: dict + nwbfile: pynwb.NWBFile + end_frame: int, optional + use for small test conversions + + """ + if nwbfile.electrodes is None: + add_electrodes_to_nwb(session_key, nwbfile) + + mapping = get_electrodes_mapping(nwbfile.electrodes) + + ecephys_module = check_module( + nwbfile, name="ecephys", description="preprocessed ephys data" + ) + + lfp = pynwb.ecephys.LFP() + ecephys_module.add(lfp) + + for ephys_recording_record in (ephys.EphysRecording & session_key).fetch( + as_dict=True + ): + probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") + + # relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") + # file_path = find_full_path(get_ephys_root_data_dir(), relative_path) + file_path = "../inbox/subject5/session1/probe_1/npx_g0_t0.imec.lf.bin" + + if ephys_recording_record["acq_software"] == "SpikeGLX": + extractor = SpikeGLXRecordingExtractor(file_path=file_path) + else: + raise ValueError( + f"unsupported acq_software type: {ephys_recording_record['acq_software']}" + ) + + if end_frame is not None: + extractor = SubRecordingExtractor(extractor, end_frame=end_frame) + + recording_channels_by_id = ( + probe.ElectrodeConfig.Electrode() & ephys_recording_record + ).fetch("electrode") + + if all(extractor.get_channel_gains() == 1): + channel_conversion = None + else: + channel_conversion = extractor.get_channel_gains() + + lfp.add_electrical_series( + pynwb.ecephys.ElectricalSeries( + name=f"ElectricalSeries{ephys_recording_record['insertion_number']}", + description=str(ephys_recording_record), + data=H5DataIO( + RecordingExtractorDataChunkIterator(extractor), compression=True + ), + rate=extractor.get_sampling_frequency(), + starting_time=( + ephys_recording_record["recording_datetime"] + - ephys_recording_record["session_datetime"] + ).total_seconds(), + electrodes=nwbfile.create_electrode_table_region( + region=[mapping[(probe_id, x)] for x in recording_channels_by_id], + name="electrodes", + description="recorded electrodes", + ), + conversion=1e-6, + channel_conversion=channel_conversion, + ) + ) + + +def ecephys_to_nwb( + session_key, + subject_id=None, + raw=True, + spikes=True, + lfp="source", + end_frame=None, +): + """ + Main function for converting ephys data to NWB + + Parameters + ---------- + session_key: dict + subject_id: str + subject_id, used if it cannot be automatically inferred + raw: bool + Whether to include the raw data from source. SpikeGLX and OpenEphys are supported + spikes: bool + Whether to include CuratedClustering + lfp: + "dj" - read LFP data from ephys.LFP + "source" - read LFP data from source (SpikeGLX supported) + False - do not convert LFP + end_frame: int, optional + Used to create small test conversions where large datasets are truncated. + """ nwbfile = session.export.nwb.session_to_nwb(session_key, subject_id=subject_id) @@ -370,15 +516,30 @@ def session_to_nwb(session_key, subject_id=None, raw=True, spikes=True, lfp=True if spikes: add_ephys_units_to_nwb(session_key, nwbfile) - if lfp: - add_ephys_lfp_to_nwb(session_key, nwbfile) + if lfp == "dj": + add_ephys_lfp_from_dj_to_nwb(session_key, nwbfile) + + if lfp == "source": + add_ephys_lfp_from_source_to_nwb(session_key, nwbfile, end_frame=end_frame) return nwbfile + def write_nwb(nwbfile, fname, check_read=True): - with NWBHDF5IO(fname, 'w') as io: + """ + Export NWBFile + + Parameters + ---------- + nwbfile: pynwb.NWBFile + fname: str + check_read: bool + If True, PyNWB will try to read the produced NWB file and ensure that it can be + read. + """ + with pynwb.NWBHDF5IO(fname, "w") as io: io.write(nwbfile) if check_read: - with NWBHDF5IO(fname, 'r') as io: + with pynwb.NWBHDF5IO(fname, "r") as io: io.read() From 470b20aa3914fb8c0bed943c7d2f82c0df5dd4fd Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Wed, 12 Jan 2022 15:41:34 -0500 Subject: [PATCH 04/77] Update element_array_ephys/export/nwb.py Co-authored-by: Dimitri Yatsenko --- element_array_ephys/export/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index dd85e8c8..1c7b321a 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -113,7 +113,7 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): description=this_probe.get("probe_comment", None), manufacturer=this_probe.get("probe_type", None), ) - shank_ids = np.unique((probe.ProbeType.Electrode & this_probe).fetch("shank")) + shank_ids = set((probe.ProbeType.Electrode & this_probe).fetch("shank")) for shank_id in shank_ids: electrode_group = nwbfile.create_electrode_group( name=f"probe{this_probe['probe']}_shank{shank_id}", From 3298341553063a8e2d165c25cfa18c0050047cde Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Wed, 12 Jan 2022 15:42:10 -0500 Subject: [PATCH 05/77] Update element_array_ephys/export/nwb.py Co-authored-by: Dimitri Yatsenko --- element_array_ephys/export/nwb.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index 1c7b321a..d5a54074 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -103,8 +103,9 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): ): insertion_record = (ephys.InsertionLocation & this_probe).fetch1() if insertion_record: - [insertion_record.pop(k) for k in ephys.InsertionLocation.heading.primary_key] - insert_location = json.dumps(insertion_record, cls=DecimalEncoder) + insert_location = json.dumps( + {k, v: for insertion_record.items() + if k not in ephys.InsertionLocation.primary_key}, cls=DecimalEncoder) else: insert_location = "unknown" From e83745247b09b369f6992590446c7167c32bddb7 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Wed, 12 Jan 2022 15:42:24 -0500 Subject: [PATCH 06/77] Update element_array_ephys/export/nwb.py Co-authored-by: Dimitri Yatsenko --- element_array_ephys/export/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index d5a54074..d6341384 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -175,7 +175,7 @@ def create_units_table( units_table.add_column( name=units_query.heading.attributes[additional_attribute].name, description=units_query.heading.attributes[additional_attribute].comment, - index=True if additional_attribute == "spike_depths" else False, + index=additional_attribute == "spike_depths", ) clustering_query = ( From d2b93f2054c51848a20d1a8da44e3cc7a9586973 Mon Sep 17 00:00:00 2001 From: bendichter Date: Fri, 14 Jan 2022 15:50:50 -0500 Subject: [PATCH 07/77] add documentation --- element_array_ephys/export/nwb.py | 44 ++++++++++++++++++++----------- 1 file changed, 29 insertions(+), 15 deletions(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index dd85e8c8..d349b67c 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -74,6 +74,19 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): Add electrodes table to NWBFile. This is needed for any ElectricalSeries, including raw source data and LFP. + ephys.InsertionLocation -> ElectrodeGroup.location + + probe.Probe::probe -> device.name + probe.Probe::probe_comment -> device.description + probe.Probe::probe_type -> device.manufacturer + + probe.ProbeType.Electrode::electrode -> electrodes["id_in_probe"] + probe.ProbeType.Electrode::y_coord -> electrodes["rel_y"] + probe.ProbeType.Electrode::x_coord -> electrodes["rel_x"] + probe.ProbeType.Electrode::shank -> electrodes["shank"] + probe.ProbeType.Electrode::shank_col -> electrodes["shank_col"] + probe.ProbeType.Electrode::shank_row -> electrodes["shank_row"] + Parameters ---------- session_key: dict @@ -82,10 +95,6 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): electrodes_query = probe.ProbeType.Electrode * probe.ElectrodeConfig.Electrode - session_electrodes = ( - probe.ElectrodeConfig.Electrode & (ephys.EphysRecording & session_key) - ).fetch() - for additional_attribute in ["shank_col", "shank_row", "shank"]: nwbfile.add_electrode_column( name=electrodes_query.heading.attributes[additional_attribute].name, @@ -148,18 +157,25 @@ def create_units_table( session_key: dict, nwbfile: pynwb.NWBFile, units_query: ephys.CuratedClustering.Unit, - paramset_record, + paramset_record: dict, name="units", desc="data on spiking units", ): """ + ephys.CuratedClustering.Unit::unit -> units.id + ephys.CuratedClustering.Unit::spike_times -> units["spike_times"] + ephys.CuratedClustering.Unit::spike_depths -> units["spike_depths"] + ephys.CuratedClustering.Unit::cluster_quality_label -> units["cluster_quality_label"] + + ephys.WaveformSet.PeakWaveform::peak_electrode_waveform -> units["waveform_mean"] + Parameters ---------- session_key: dict nwbfile: pynwb.NWBFile units_query: ephys.CuratedClustering.Unit - paramset_record: int + paramset_record: dict name: str, optional default="units" desc: str, optional @@ -223,8 +239,6 @@ def add_ephys_units_to_nwb( session_key: dict nwbfile: pynwb.NWBFile primary_clustering_paramset_idx: int, optional - - """ if not ephys.ClusteringTask & session_key: @@ -252,6 +266,7 @@ def add_ephys_units_to_nwb( name = f"units_{paramset_record['clustering_method']}" units_table = create_units_table( session_key, + nwbfile, units_query, paramset_record, name=name, @@ -295,9 +310,6 @@ def add_ephys_recording_to_nwb( Used for small test conversions """ - if not ephys.EphysRecording & session_key: - return - if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) @@ -360,15 +372,15 @@ def add_ephys_lfp_from_dj_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): """ Read LFP data from the data in element-aray-ephys + ephys.LFP::lfp -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data + ephys.LFP::lfp_time_stamps -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].timestamps + Parameters ---------- session_key: dict nwbfile: NWBFile """ - if not ephys.LFP & session_key: - return - if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) @@ -407,7 +419,9 @@ def add_ephys_lfp_from_source_to_nwb( session_key: dict, nwbfile: pynwb.NWBFile, end_frame=None ): """ - Read the LFP data directly from the source file. Currently only works for SpikeGLX data. + Read the LFP data directly from the source file. Currently, only works for SpikeGLX data. + + ephys.EphysRecording::recording_datetime -> acquisition Parameters ---------- From d924d57ea62ec24c791135a16f3ed397493084d5 Mon Sep 17 00:00:00 2001 From: bendichter Date: Fri, 14 Jan 2022 16:02:36 -0500 Subject: [PATCH 08/77] * optimize imports * black * upgrade to latest version of conversion-tools * upgrade to latest spikeinterface api * --- element_array_ephys/export/nwb.py | 145 ++++++++++++++++-------------- 1 file changed, 77 insertions(+), 68 deletions(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index d349b67c..ef239a39 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -1,21 +1,18 @@ import decimal +import importlib +import json import numpy as np +import pynwb +from element_session import session from hdmf.backends.hdf5 import H5DataIO +from nwb_conversion_tools.utils.conversion_tools import get_module from nwb_conversion_tools.utils.genericdatachunkiterator import GenericDataChunkIterator -from nwb_conversion_tools.utils.recordingextractordatachunkiterator import ( - RecordingExtractorDataChunkIterator, -) -from nwb_conversion_tools.utils.spike_interface import check_module - -import pynwb -from spikeextractors import ( - SpikeGLXRecordingExtractor, - SubRecordingExtractor, - OpenEphysNPIXRecordingExtractor, +from nwb_conversion_tools.utils.spikeinterfacerecordingdatachunkiterator import ( + SpikeInterfaceRecordingDataChunkIterator, ) +from spikeinterface import extractors from tqdm import tqdm - -from element_session import session +from uuid import uuid4 # from element_data_loader.utils import find_full_path @@ -46,7 +43,7 @@ def __init__(self, lfp_electrodes_query, chunk_length: int = 10000): self.electrodes = self.lfp_electrodes_query.fetch("electrode") first_record = ( - self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) + self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) ).fetch(as_dict=True)[0] self.n_channels = len(lfp_electrodes_query) @@ -108,11 +105,14 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): ) for this_probe in (ephys.ProbeInsertion * probe.Probe & session_key).fetch( - as_dict=True + as_dict=True ): insertion_record = (ephys.InsertionLocation & this_probe).fetch1() if insertion_record: - [insertion_record.pop(k) for k in ephys.InsertionLocation.heading.primary_key] + [ + insertion_record.pop(k) + for k in ephys.InsertionLocation.heading.primary_key + ] insert_location = json.dumps(insertion_record, cls=DecimalEncoder) else: insert_location = "unknown" @@ -132,7 +132,7 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): ) electrodes_query = ( - probe.ProbeType.Electrode & this_probe & dict(shank=shank_id) + probe.ProbeType.Electrode & this_probe & dict(shank=shank_id) ).fetch(as_dict=True) for electrode in electrodes_query: nwbfile.add_electrode( @@ -154,12 +154,12 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): def create_units_table( - session_key: dict, - nwbfile: pynwb.NWBFile, - units_query: ephys.CuratedClustering.Unit, - paramset_record: dict, - name="units", - desc="data on spiking units", + session_key: dict, + nwbfile: pynwb.NWBFile, + units_query: ephys.CuratedClustering.Unit, + paramset_record, + name="units", + desc="data on spiking units", ): """ @@ -175,7 +175,7 @@ def create_units_table( session_key: dict nwbfile: pynwb.NWBFile units_query: ephys.CuratedClustering.Unit - paramset_record: dict + paramset_record: int name: str, optional default="units" desc: str, optional @@ -194,23 +194,23 @@ def create_units_table( ) clustering_query = ( - ephys.EphysRecording * ephys.ClusteringTask & session_key & paramset_record + ephys.EphysRecording * ephys.ClusteringTask & session_key & paramset_record ) for unit in tqdm( - (clustering_query @ ephys.CuratedClustering.Unit).fetch(as_dict=True), - desc=f"creating units table for paramset {paramset_record['paramset_idx']}", + (clustering_query @ ephys.CuratedClustering.Unit).fetch(as_dict=True), + desc=f"creating units table for paramset {paramset_record['paramset_idx']}", ): probe_id, shank_num = ( - probe.ProbeType.Electrode - * ephys.ProbeInsertion - * ephys.CuratedClustering.Unit - & {"unit": unit["unit"], "insertion_number": unit["insertion_number"]} + probe.ProbeType.Electrode + * ephys.ProbeInsertion + * ephys.CuratedClustering.Unit + & {"unit": unit["unit"], "insertion_number": unit["insertion_number"]} ).fetch1("probe", "shank") waveform_mean = ( - ephys.WaveformSet.PeakWaveform() & clustering_query & unit + ephys.WaveformSet.PeakWaveform() & clustering_query & unit ).fetch1("peak_electrode_waveform") units_table.add_row( @@ -229,7 +229,7 @@ def create_units_table( def add_ephys_units_to_nwb( - session_key: dict, nwbfile: pynwb.NWBFile, primary_clustering_paramset_idx: int = 0 + session_key: dict, nwbfile: pynwb.NWBFile, primary_clustering_paramset_idx: int = 0 ): """ Add spiking data to NWBFile. @@ -251,7 +251,7 @@ def add_ephys_units_to_nwb( units_query = ephys.CuratedClustering.Unit() & session_key for paramset_record in ( - ephys.ClusteringParamSet & ephys.CuratedClustering & session_key + ephys.ClusteringParamSet & ephys.CuratedClustering & session_key ).fetch("paramset_idx", "clustering_method", "paramset_desc", as_dict=True): if paramset_record["paramset_idx"] == primary_clustering_paramset_idx: units_table = create_units_table( @@ -266,13 +266,12 @@ def add_ephys_units_to_nwb( name = f"units_{paramset_record['clustering_method']}" units_table = create_units_table( session_key, - nwbfile, units_query, paramset_record, name=name, desc=paramset_record["paramset_desc"], ) - ecephys_module = check_module(nwbfile, "ecephys") + ecephys_module = get_module(nwbfile, "ecephys") ecephys_module.add(units_table) @@ -298,7 +297,7 @@ def get_electrodes_mapping(electrodes): def add_ephys_recording_to_nwb( - session_key: dict, nwbfile: pynwb.NWBFile, end_frame: int = None + session_key: dict, nwbfile: pynwb.NWBFile, end_frame: int = None ): """ @@ -316,18 +315,22 @@ def add_ephys_recording_to_nwb( mapping = get_electrodes_mapping(nwbfile.electrodes) for ephys_recording_record in (ephys.EphysRecording & session_key).fetch( - as_dict=True + as_dict=True ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") # relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") # file_path = find_full_path(get_ephys_root_data_dir(), relative_path) - file_path = "../inbox/subject5/session1/probe_1/npx_g0_t0.imec.ap.bin" + + # ephys_recording_record["acq_software"] = "OpenEphys" + # file_path = "../sciopsdemo-inbox/allen/SixProbesExperiments" + + file_path = "../inbox/subject5/session1/probe_1/" if ephys_recording_record["acq_software"] == "SpikeGLX": - extractor = SpikeGLXRecordingExtractor(file_path=file_path) + extractor = extractors.read_spikeglx(file_path, "imec.ap") elif ephys_recording_record["acq_software"] == "OpenEphys": - extractor = OpenEphysNPIXRecordingExtractor(file_path=file_path) + extractor = extractors.read_openephys(file_path, stream_id="0") else: raise ValueError( f"unsupported acq_software type: {ephys_recording_record['acq_software']}" @@ -339,23 +342,21 @@ def add_ephys_recording_to_nwb( channel_conversion = extractor.get_channel_gains() if end_frame is not None: - extractor = SubRecordingExtractor(extractor, end_frame=end_frame) + extractor = extractor.frame_slice(0, end_frame) recording_channels_by_id = ( - probe.ElectrodeConfig.Electrode() & ephys_recording_record + probe.ElectrodeConfig.Electrode() & ephys_recording_record ).fetch("electrode") nwbfile.add_acquisition( pynwb.ecephys.ElectricalSeries( name=f"ElectricalSeries{ephys_recording_record['insertion_number']}", description=str(ephys_recording_record), - data=H5DataIO( - RecordingExtractorDataChunkIterator(extractor), compression=True - ), + data=SpikeInterfaceRecordingDataChunkIterator(extractor), rate=ephys_recording_record["sampling_rate"], starting_time=( - ephys_recording_record["recording_datetime"] - - ephys_recording_record["session_datetime"] + ephys_recording_record["recording_datetime"] + - ephys_recording_record["session_datetime"] ).total_seconds(), electrodes=nwbfile.create_electrode_table_region( region=[mapping[(probe_id, x)] for x in recording_channels_by_id], @@ -384,7 +385,7 @@ def add_ephys_lfp_from_dj_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) - ecephys_module = check_module( + ecephys_module = get_module( nwbfile, name="ecephys", description="preprocessed ephys data" ) @@ -416,7 +417,7 @@ def add_ephys_lfp_from_dj_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): def add_ephys_lfp_from_source_to_nwb( - session_key: dict, nwbfile: pynwb.NWBFile, end_frame=None + session_key: dict, nwbfile: pynwb.NWBFile, end_frame=None ): """ Read the LFP data directly from the source file. Currently, only works for SpikeGLX data. @@ -436,7 +437,7 @@ def add_ephys_lfp_from_source_to_nwb( mapping = get_electrodes_mapping(nwbfile.electrodes) - ecephys_module = check_module( + ecephys_module = get_module( nwbfile, name="ecephys", description="preprocessed ephys data" ) @@ -444,26 +445,26 @@ def add_ephys_lfp_from_source_to_nwb( ecephys_module.add(lfp) for ephys_recording_record in (ephys.EphysRecording & session_key).fetch( - as_dict=True + as_dict=True ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") # relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") # file_path = find_full_path(get_ephys_root_data_dir(), relative_path) - file_path = "../inbox/subject5/session1/probe_1/npx_g0_t0.imec.lf.bin" + file_path = "../inbox/subject5/session1/probe_1/" if ephys_recording_record["acq_software"] == "SpikeGLX": - extractor = SpikeGLXRecordingExtractor(file_path=file_path) + extractor = extractors.read_spikeglx(file_path, "imec.lf") else: raise ValueError( f"unsupported acq_software type: {ephys_recording_record['acq_software']}" ) if end_frame is not None: - extractor = SubRecordingExtractor(extractor, end_frame=end_frame) + extractor = extractor.frame_slice(0, end_frame) recording_channels_by_id = ( - probe.ElectrodeConfig.Electrode() & ephys_recording_record + probe.ElectrodeConfig.Electrode() & ephys_recording_record ).fetch("electrode") if all(extractor.get_channel_gains() == 1): @@ -475,13 +476,11 @@ def add_ephys_lfp_from_source_to_nwb( pynwb.ecephys.ElectricalSeries( name=f"ElectricalSeries{ephys_recording_record['insertion_number']}", description=str(ephys_recording_record), - data=H5DataIO( - RecordingExtractorDataChunkIterator(extractor), compression=True - ), + data=SpikeInterfaceRecordingDataChunkIterator(extractor), rate=extractor.get_sampling_frequency(), starting_time=( - ephys_recording_record["recording_datetime"] - - ephys_recording_record["session_datetime"] + ephys_recording_record["recording_datetime"] + - ephys_recording_record["session_datetime"] ).total_seconds(), electrodes=nwbfile.create_electrode_table_region( region=[mapping[(probe_id, x)] for x in recording_channels_by_id], @@ -494,13 +493,15 @@ def add_ephys_lfp_from_source_to_nwb( ) -def ecephys_to_nwb( - session_key, - subject_id=None, - raw=True, - spikes=True, - lfp="source", - end_frame=None, +def ecephys_session_to_nwb( + session_key, + session_identifier=str(uuid4()), + session_start_time=datetime.datetime.now(), + subject_id=None, + raw=True, + spikes=True, + lfp="source", + end_frame=None, ): """ Main function for converting ephys data to NWB @@ -508,6 +509,8 @@ def ecephys_to_nwb( Parameters ---------- session_key: dict + session_identifier: str, optional + session_start_time: datetime.datetime, optional subject_id: str subject_id, used if it cannot be automatically inferred raw: bool @@ -522,7 +525,13 @@ def ecephys_to_nwb( Used to create small test conversions where large datasets are truncated. """ - nwbfile = session.export.nwb.session_to_nwb(session_key, subject_id=subject_id) + if session.schema.is_activated(): + nwbfile = session.export.nwb.session_to_nwb(session_key, subject_id=subject_id) + else: + nwbfile = pynwb.NWBFile( + "automatically generated", session_identifier, session_start_time + ) + nwbfile.subject = pynwb.file.Subject(subject_id=subject_id) if raw: add_ephys_recording_to_nwb(session_key, nwbfile, end_frame=end_frame) From a7f846be82b020c5b6ace34b524dd99c5438d345 Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 25 Jan 2022 14:36:51 -0500 Subject: [PATCH 09/77] add missing nwbfile arg --- element_array_ephys/export/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index c86bed0d..51a5f3a7 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -265,6 +265,7 @@ def add_ephys_units_to_nwb( name = f"units_{paramset_record['clustering_method']}" units_table = create_units_table( session_key, + nwbfile, units_query, paramset_record, name=name, From b35ee723dc39777dfede2334452021a5c234f6d8 Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 25 Jan 2022 15:01:07 -0500 Subject: [PATCH 10/77] add some docstrings --- element_array_ephys/export/nwb.py | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb.py index 51a5f3a7..0ecde99c 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb.py @@ -233,6 +233,13 @@ def add_ephys_units_to_nwb( """ Add spiking data to NWBFile. + ephys.CuratedClustering.Unit::unit -> units.id + ephys.CuratedClustering.Unit::spike_times -> units["spike_times"] + ephys.CuratedClustering.Unit::spike_depths -> units["spike_depths"] + ephys.CuratedClustering.Unit::cluster_quality_label -> units["cluster_quality_label"] + + ephys.WaveformSet.PeakWaveform::peak_electrode_waveform -> units["waveform_mean"] + Parameters ---------- session_key: dict @@ -299,7 +306,11 @@ def get_electrodes_mapping(electrodes): def add_ephys_recording_to_nwb( session_key: dict, nwbfile: pynwb.NWBFile, end_frame: int = None ): - """ + """Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically + applies lossless compression to the data, to the final file might be smaller than the original, but there is no + data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. + + source data -> acquisition["ElectricalSeries"] Parameters ---------- From 19f78b78ef4467f1edd7b3270a8df51cdf349a91 Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 25 Jan 2022 15:01:28 -0500 Subject: [PATCH 11/77] add readme for exporting to NWB --- element_array_ephys/export/export_nwb.md | 55 ++++++++++++++++++++++++ 1 file changed, 55 insertions(+) create mode 100644 element_array_ephys/export/export_nwb.md diff --git a/element_array_ephys/export/export_nwb.md b/element_array_ephys/export/export_nwb.md new file mode 100644 index 00000000..44e90fd3 --- /dev/null +++ b/element_array_ephys/export/export_nwb.md @@ -0,0 +1,55 @@ +# Exporting data to NWB + +The export/nwb.py module maps from the element-array-ephys data structure to NWB. +The main function is `ecephys_session_to_nwb`, which contains flags to control calling the following functions, +which can be called independently: + +1. `session.export.nwb.session_to_nwb`: Gathers session-level metadata + + +2. `add_electrodes_to_nwb`: Add electrodes table to NWBFile. This is needed for any ElectricalSeries, including + raw source data and LFP. + + + ephys.InsertionLocation -> ElectrodeGroup.location + + probe.Probe::probe -> device.name + probe.Probe::probe_comment -> device.description + probe.Probe::probe_type -> device.manufacturer + + probe.ProbeType.Electrode::electrode -> electrodes["id_in_probe"] + probe.ProbeType.Electrode::y_coord -> electrodes["rel_y"] + probe.ProbeType.Electrode::x_coord -> electrodes["rel_x"] + probe.ProbeType.Electrode::shank -> electrodes["shank"] + probe.ProbeType.Electrode::shank_col -> electrodes["shank_col"] + probe.ProbeType.Electrode::shank_row -> electrodes["shank_row"] + +3. `add_ephys_recording_to_nwb`: Read voltage data directly from source files and iteratively transfer them to the + NWB file. Automatically applies lossless compression to the data, to the final file might be smaller than the original, but there is no + data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. + + + source data -> acquisition["ElectricalSeries"] + +4. `add_ephys_units_to_nwb`: Add spiking data from CuratedClustering to NWBFile. + + + ephys.CuratedClustering.Unit::unit -> units.id + ephys.CuratedClustering.Unit::spike_times -> units["spike_times"] + ephys.CuratedClustering.Unit::spike_depths -> units["spike_depths"] + ephys.CuratedClustering.Unit::cluster_quality_label -> units["cluster_quality_label"] + + ephys.WaveformSet.PeakWaveform::peak_electrode_waveform -> units["waveform_mean"] + +5. `add_ephys_lfp_from_dj_to_nwb`: Read LFP data from the data in element-array-ephys and convert to NWB. + + + ephys.LFP::lfp -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data + ephys.LFP::lfp_time_stamps -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].timestamps + +6. `add_ephys_lfp_from_source_to_nwb`: Read the LFP data directly from the source file. Currently, only works for + SpikeGLX data. Can be used instead of `add_ephys_lfp_from_dj_to_nwb`. + + + source data -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data + source data -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].timestamps \ No newline at end of file From 013ae7da1bc94ed2c6c5e87b79378e121826882b Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 15:48:00 -0500 Subject: [PATCH 12/77] refactor to include requirements for nwb conversion --- .../export/{export_nwb.md => nwb/README.md} | 2 +- element_array_ephys/export/nwb/__init__.py | 0 element_array_ephys/export/{ => nwb}/nwb.py | 23 +++++++------------ .../export/nwb/requirements.txt | 4 ++++ 4 files changed, 13 insertions(+), 16 deletions(-) rename element_array_ephys/export/{export_nwb.md => nwb/README.md} (96%) create mode 100644 element_array_ephys/export/nwb/__init__.py rename element_array_ephys/export/{ => nwb}/nwb.py (95%) create mode 100644 element_array_ephys/export/nwb/requirements.txt diff --git a/element_array_ephys/export/export_nwb.md b/element_array_ephys/export/nwb/README.md similarity index 96% rename from element_array_ephys/export/export_nwb.md rename to element_array_ephys/export/nwb/README.md index 44e90fd3..5b146785 100644 --- a/element_array_ephys/export/export_nwb.md +++ b/element_array_ephys/export/nwb/README.md @@ -1,6 +1,6 @@ # Exporting data to NWB -The export/nwb.py module maps from the element-array-ephys data structure to NWB. +The `export/nwb/nwb.py` module maps from the element-array-ephys data structure to NWB. The main function is `ecephys_session_to_nwb`, which contains flags to control calling the following functions, which can be called independently: diff --git a/element_array_ephys/export/nwb/__init__.py b/element_array_ephys/export/nwb/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/element_array_ephys/export/nwb.py b/element_array_ephys/export/nwb/nwb.py similarity index 95% rename from element_array_ephys/export/nwb.py rename to element_array_ephys/export/nwb/nwb.py index 0ecde99c..d93adba3 100644 --- a/element_array_ephys/export/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -14,8 +14,7 @@ from tqdm import tqdm from uuid import uuid4 - -# from element_data_loader.utils import find_full_path +from workflow.pipeline import ephys class DecimalEncoder(json.JSONEncoder): @@ -155,7 +154,7 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): def create_units_table( session_key: dict, nwbfile: pynwb.NWBFile, - units_query: ephys.CuratedClustering.Unit, + units_query, paramset_record, name="units", desc="data on spiking units", @@ -330,16 +329,11 @@ def add_ephys_recording_to_nwb( ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") - # relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") - # file_path = find_full_path(get_ephys_root_data_dir(), relative_path) - - # ephys_recording_record["acq_software"] = "OpenEphys" - # file_path = "../sciopsdemo-inbox/allen/SixProbesExperiments" - - file_path = "../inbox/subject5/session1/probe_1/" + relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") + file_path = ephys.find_full_path(get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": - extractor = extractors.read_spikeglx(file_path, "imec.ap") + extractor = extractors.read_spikeglx(os.path.split(file_path)[0], "imec.ap") elif ephys_recording_record["acq_software"] == "OpenEphys": extractor = extractors.read_openephys(file_path, stream_id="0") else: @@ -460,12 +454,11 @@ def add_ephys_lfp_from_source_to_nwb( ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") - # relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") - # file_path = find_full_path(get_ephys_root_data_dir(), relative_path) - file_path = "../inbox/subject5/session1/probe_1/" + relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") + file_path = find_full_path(ephys.get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": - extractor = extractors.read_spikeglx(file_path, "imec.lf") + extractor = extractors.read_spikeglx(os.path.split(file_path)[0], "imec.lf") else: raise ValueError( f"unsupported acq_software type: {ephys_recording_record['acq_software']}" diff --git a/element_array_ephys/export/nwb/requirements.txt b/element_array_ephys/export/nwb/requirements.txt new file mode 100644 index 00000000..e4c20c0b --- /dev/null +++ b/element_array_ephys/export/nwb/requirements.txt @@ -0,0 +1,4 @@ +dandi +git+https://github.com/catalystneuro/nwb-conversion-tools.git +git+https://github.com/SpikeInterface/spikeinterface.git +git+https://github.com/bendichter/python-neo.git@2d0ba5213f6e2add29dd13b2131650e3325edb23 \ No newline at end of file From 13632147a56826fe7ccba065aa02e370950a72c0 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:02:33 -0500 Subject: [PATCH 13/77] relative import of ephys --- element_array_ephys/export/nwb/nwb.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index d93adba3..aa4270ea 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -14,7 +14,8 @@ from tqdm import tqdm from uuid import uuid4 -from workflow.pipeline import ephys +from ... import ephys +#from workflow.pipeline import ephys class DecimalEncoder(json.JSONEncoder): From d1f3daba676385a54aa108bb554bce008594b96d Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:04:08 -0500 Subject: [PATCH 14/77] add datetime import --- element_array_ephys/export/nwb/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index aa4270ea..2c6bb214 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -1,3 +1,4 @@ +import datetime import decimal import importlib import json From 365b43b0fa584a531b386afe74eb55e853773a2d Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:08:06 -0500 Subject: [PATCH 15/77] add optional session keys to ecephys_session_to_nwb --- element_array_ephys/export/nwb/nwb.py | 47 ++++++++++++++++++--------- 1 file changed, 32 insertions(+), 15 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 2c6bb214..ee94d4fa 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -500,14 +500,16 @@ def add_ephys_lfp_from_source_to_nwb( def ecephys_session_to_nwb( - session_key, - session_identifier=str(uuid4()), - session_start_time=datetime.datetime.now(), - subject_id=None, - raw=True, - spikes=True, - lfp="source", - end_frame=None, + session_key, + subject_id=None, + raw=True, + spikes=True, + lfp="source", + end_frame=None, + lab_key=None, + project_key=None, + protocol_key=None, + nwbfile_kwargs=None, ): """ Main function for converting ephys data to NWB @@ -515,8 +517,6 @@ def ecephys_session_to_nwb( Parameters ---------- session_key: dict - session_identifier: str, optional - session_start_time: datetime.datetime, optional subject_id: str subject_id, used if it cannot be automatically inferred raw: bool @@ -529,15 +529,31 @@ def ecephys_session_to_nwb( False - do not convert LFP end_frame: int, optional Used to create small test conversions where large datasets are truncated. + lab_key, project_key, and protocol_key: dictionaries used to look up optional additional metadata + nwbfile_kwargs: dict, optional + - If element-session is not being used, this argument is required and must be a dictionary containing + 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime), + the minimal data for instantiating an NWBFile object. + + - If element-session is being used, this argument can optionally be used to add over overwrite NWBFile fields. """ if session.schema.is_activated(): - nwbfile = session.export.nwb.session_to_nwb(session_key, subject_id=subject_id) - else: - nwbfile = pynwb.NWBFile( - "automatically generated", session_identifier, session_start_time + nwbfile = session_to_nwb( + session_key, + subject_id=subject_id, + lab_key=lab_key, + project_key=project_key, + protocol_key=protocol_key, + additional_nwbfile_kwargs=nwbfile_kwargs, ) - nwbfile.subject = pynwb.file.Subject(subject_id=subject_id) + else: + if not isinstance(nwbfile_kwargs, dict) and {'session_description', 'identifier', 'session_start_time'}.issubset(nwbfile_kwargs): + raise ValueError( + "If element-session is not activated, you must include nwbfile_kwargs as a dictionary." + "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)" + ) + nwbfile = pynwb.NWBFile(**nwbfile_kwargs) if raw: add_ephys_recording_to_nwb(session_key, nwbfile, end_frame=end_frame) @@ -554,6 +570,7 @@ def ecephys_session_to_nwb( return nwbfile + def write_nwb(nwbfile, fname, check_read=True): """ Export NWBFile From 96c57f756243eb00ced7e157e417ab1af5881c10 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:10:16 -0500 Subject: [PATCH 16/77] fix nwbfile_kwargs logic --- element_array_ephys/export/nwb/nwb.py | 30 +++++++++++++++------------ 1 file changed, 17 insertions(+), 13 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index ee94d4fa..8c1b12d5 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -16,7 +16,8 @@ from uuid import uuid4 from ... import ephys -#from workflow.pipeline import ephys + +# from workflow.pipeline import ephys class DecimalEncoder(json.JSONEncoder): @@ -500,16 +501,16 @@ def add_ephys_lfp_from_source_to_nwb( def ecephys_session_to_nwb( - session_key, - subject_id=None, - raw=True, - spikes=True, - lfp="source", - end_frame=None, - lab_key=None, - project_key=None, - protocol_key=None, - nwbfile_kwargs=None, + session_key, + subject_id=None, + raw=True, + spikes=True, + lfp="source", + end_frame=None, + lab_key=None, + project_key=None, + protocol_key=None, + nwbfile_kwargs=None, ): """ Main function for converting ephys data to NWB @@ -548,7 +549,11 @@ def ecephys_session_to_nwb( additional_nwbfile_kwargs=nwbfile_kwargs, ) else: - if not isinstance(nwbfile_kwargs, dict) and {'session_description', 'identifier', 'session_start_time'}.issubset(nwbfile_kwargs): + if isinstance(nwbfile_kwargs, dict) and not { + "session_description", + "identifier", + "session_start_time", + }.issubset(nwbfile_kwargs): raise ValueError( "If element-session is not activated, you must include nwbfile_kwargs as a dictionary." "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)" @@ -570,7 +575,6 @@ def ecephys_session_to_nwb( return nwbfile - def write_nwb(nwbfile, fname, check_read=True): """ Export NWBFile From eb47ee506db1fd1929da2d557c17776c28675326 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:12:38 -0500 Subject: [PATCH 17/77] fix nwbfile_kwargs logic --- element_array_ephys/export/nwb/nwb.py | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 8c1b12d5..4694bcaf 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -549,11 +549,10 @@ def ecephys_session_to_nwb( additional_nwbfile_kwargs=nwbfile_kwargs, ) else: - if isinstance(nwbfile_kwargs, dict) and not { - "session_description", - "identifier", - "session_start_time", - }.issubset(nwbfile_kwargs): + if not ( + isinstance(nwbfile_kwargs, dict) + and {"session_description", "identifier", "session_start_time"}.issubset(nwbfile_kwargs) + ): raise ValueError( "If element-session is not activated, you must include nwbfile_kwargs as a dictionary." "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)" From 1b676293227e03e77293c223f8d3f336c832ee59 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:14:46 -0500 Subject: [PATCH 18/77] import from workflow pipeline --- element_array_ephys/export/nwb/nwb.py | 20 ++++++++++++++------ 1 file changed, 14 insertions(+), 6 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 4694bcaf..a927dc2b 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -15,9 +15,9 @@ from tqdm import tqdm from uuid import uuid4 -from ... import ephys +#from ... import ephys -# from workflow.pipeline import ephys +from workflow.pipeline import ephys class DecimalEncoder(json.JSONEncoder): @@ -112,8 +112,12 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): insertion_record = (ephys.InsertionLocation & this_probe).fetch1() if insertion_record: insert_location = json.dumps( - {k: v for k, v in insertion_record.items() if k not in ephys.InsertionLocation.primary_key}, - cls=DecimalEncoder + { + k: v + for k, v in insertion_record.items() + if k not in ephys.InsertionLocation.primary_key + }, + cls=DecimalEncoder, ) else: insert_location = "unknown" @@ -332,7 +336,9 @@ def add_ephys_recording_to_nwb( ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") - relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") + relative_path = ( + ephys.EphysRecording.EphysFile & ephys_recording_record + ).fetch1("file_path") file_path = ephys.find_full_path(get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": @@ -457,7 +463,9 @@ def add_ephys_lfp_from_source_to_nwb( ): probe_id = (ephys.ProbeInsertion() & ephys_recording_record).fetch1("probe") - relative_path = (ephys.EphysRecording.EphysFile & ephys_recording_record).fetch1("file_path") + relative_path = ( + ephys.EphysRecording.EphysFile & ephys_recording_record + ).fetch1("file_path") file_path = find_full_path(ephys.get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": From 5de22e06fa2db7b6e2005bbc094b05fc3581e850 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:17:14 -0500 Subject: [PATCH 19/77] import from workflow pipeline --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index a927dc2b..0a008b4d 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -17,7 +17,7 @@ #from ... import ephys -from workflow.pipeline import ephys +from workflow_array_ephys.pipeline import ephys class DecimalEncoder(json.JSONEncoder): From 1dccb6893fa8a77408e598f278995998b5c8b45c Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:18:56 -0500 Subject: [PATCH 20/77] import session_to_nwb --- element_array_ephys/export/nwb/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 0a008b4d..69625436 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -548,6 +548,7 @@ def ecephys_session_to_nwb( """ if session.schema.is_activated(): + from element_session.export.nwb import session_to_nwb nwbfile = session_to_nwb( session_key, subject_id=subject_id, From 379ae11c718242d04cc6c2eb62a76a433d2b94cf Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:40:10 -0500 Subject: [PATCH 21/77] import probe --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 69625436..ac26054b 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -17,7 +17,7 @@ #from ... import ephys -from workflow_array_ephys.pipeline import ephys +from workflow_array_ephys.pipeline import ephys, probe class DecimalEncoder(json.JSONEncoder): From e99247874924a39596510fb144f3272f00886804 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:42:52 -0500 Subject: [PATCH 22/77] ephys.get_ephys_root_data_dir --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index ac26054b..c184a448 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -339,7 +339,7 @@ def add_ephys_recording_to_nwb( relative_path = ( ephys.EphysRecording.EphysFile & ephys_recording_record ).fetch1("file_path") - file_path = ephys.find_full_path(get_ephys_root_data_dir(), relative_path) + file_path = ephys.find_full_path(ephys.get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": extractor = extractors.read_spikeglx(os.path.split(file_path)[0], "imec.ap") From 8d3df711fc6f79282e81ac2794c2ca0b5d19c10c Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:51:18 -0500 Subject: [PATCH 23/77] standardize slashes --- element_array_ephys/export/nwb/nwb.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index c184a448..e0dac593 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -339,6 +339,7 @@ def add_ephys_recording_to_nwb( relative_path = ( ephys.EphysRecording.EphysFile & ephys_recording_record ).fetch1("file_path") + relative_path = relative_path.replace("\\","/") file_path = ephys.find_full_path(ephys.get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": @@ -466,6 +467,7 @@ def add_ephys_lfp_from_source_to_nwb( relative_path = ( ephys.EphysRecording.EphysFile & ephys_recording_record ).fetch1("file_path") + relative_path = relative_path.replace("\\","/") file_path = find_full_path(ephys.get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": From ba3f86a2a956e26beb7223cbf09dc779de894a8b Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:52:02 -0500 Subject: [PATCH 24/77] import os --- element_array_ephys/export/nwb/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index e0dac593..b8c6ccdf 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -1,3 +1,4 @@ +import os import datetime import decimal import importlib From 654d567f8701940bdec6c426d65cb619c3fa2016 Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 30 Jan 2022 17:53:14 -0500 Subject: [PATCH 25/77] ephys.find_full_path --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index b8c6ccdf..80778090 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -469,7 +469,7 @@ def add_ephys_lfp_from_source_to_nwb( ephys.EphysRecording.EphysFile & ephys_recording_record ).fetch1("file_path") relative_path = relative_path.replace("\\","/") - file_path = find_full_path(ephys.get_ephys_root_data_dir(), relative_path) + file_path = ephys.find_full_path(ephys.get_ephys_root_data_dir(), relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": extractor = extractors.read_spikeglx(os.path.split(file_path)[0], "imec.lf") From 7deb00f8e360c5b23e8b63333c55b50ece5334c3 Mon Sep 17 00:00:00 2001 From: bendichter Date: Fri, 4 Feb 2022 14:58:03 -0400 Subject: [PATCH 26/77] correctly set conversion and channel_conversion --- element_array_ephys/export/nwb/nwb.py | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 80778090..785f4d7e 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -485,10 +485,14 @@ def add_ephys_lfp_from_source_to_nwb( probe.ElectrodeConfig.Electrode() & ephys_recording_record ).fetch("electrode") - if all(extractor.get_channel_gains() == 1): - channel_conversion = None - else: - channel_conversion = extractor.get_channel_gains() + conversion = 1e-6 + channel_conversion = None + if any(extractor.get_channel_gains() != 1): + gains = extractor.get_channel_gains() + if all(x == gains[0] for x in gains): + conversion *= gains[0] + else: + channel_conversion = gains lfp.add_electrical_series( pynwb.ecephys.ElectricalSeries( @@ -505,7 +509,7 @@ def add_ephys_lfp_from_source_to_nwb( name="electrodes", description="recorded electrodes", ), - conversion=1e-6, + conversion=conversion, channel_conversion=channel_conversion, ) ) From 79536621e8b6b9f35ed4cb61e5384c433da28bc0 Mon Sep 17 00:00:00 2001 From: bendichter Date: Fri, 4 Feb 2022 15:09:06 -0400 Subject: [PATCH 27/77] refactor into gains_helper --- element_array_ephys/export/nwb/nwb.py | 28 ++++++++++++--------------- 1 file changed, 12 insertions(+), 16 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 785f4d7e..372ad02c 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -310,6 +310,14 @@ def get_electrodes_mapping(electrodes): } +def gains_helper(gains): + if all(x == 1 for x in gains): + return dict(conversion=1e-6, channel_conversion=None) + if all(x == gains[0] for x in gains): + return dict(conversion=1e-6 * gains[0], channel_conversion=None) + return dict(conversion=1e-6, channel_conversion=gains) + + def add_ephys_recording_to_nwb( session_key: dict, nwbfile: pynwb.NWBFile, end_frame: int = None ): @@ -352,10 +360,7 @@ def add_ephys_recording_to_nwb( f"unsupported acq_software type: {ephys_recording_record['acq_software']}" ) - if all(extractor.get_channel_gains() == 1): - channel_conversion = None - else: - channel_conversion = extractor.get_channel_gains() + conversion_kwargs = gains_helper(extractor.get_channel_gains()) if end_frame is not None: extractor = extractor.frame_slice(0, end_frame) @@ -379,8 +384,7 @@ def add_ephys_recording_to_nwb( name="electrodes", description="recorded electrodes", ), - conversion=1e-6, - channel_conversion=channel_conversion, + **conversion_kwargs ) ) @@ -485,14 +489,7 @@ def add_ephys_lfp_from_source_to_nwb( probe.ElectrodeConfig.Electrode() & ephys_recording_record ).fetch("electrode") - conversion = 1e-6 - channel_conversion = None - if any(extractor.get_channel_gains() != 1): - gains = extractor.get_channel_gains() - if all(x == gains[0] for x in gains): - conversion *= gains[0] - else: - channel_conversion = gains + conversion_kwargs = gains_helper(extractor.get_channel_gains()) lfp.add_electrical_series( pynwb.ecephys.ElectricalSeries( @@ -509,8 +506,7 @@ def add_ephys_lfp_from_source_to_nwb( name="electrodes", description="recorded electrodes", ), - conversion=conversion, - channel_conversion=channel_conversion, + **conversion_kwargs ) ) From 4cab8c86620094e912e7193c3a579989046d26b1 Mon Sep 17 00:00:00 2001 From: bendichter Date: Fri, 4 Feb 2022 16:00:20 -0400 Subject: [PATCH 28/77] add tests --- element_array_ephys/export/nwb/tests.py | 37 +++++++++++++++++++++++++ 1 file changed, 37 insertions(+) create mode 100644 element_array_ephys/export/nwb/tests.py diff --git a/element_array_ephys/export/nwb/tests.py b/element_array_ephys/export/nwb/tests.py new file mode 100644 index 00000000..4525038a --- /dev/null +++ b/element_array_ephys/export/nwb/tests.py @@ -0,0 +1,37 @@ +from element_array_ephys.export.nwb.nwb import ecephys_session_to_nwb, write_nwb + +from pynwb.ecephys import ElectricalSeries + + +def test_convert_to_nwb(): + + nwbfile = ecephys_session_to_nwb( + dict(subject="subject5", session_datetime="2020-05-12 04:13:07") + ) + + for x in ("262716621", "714000838"): + assert x in nwbfile.devices + + assert len(nwbfile.electrodes) == 1920 + for col in ("shank", "shank_row", "shank_col"): + assert col in nwbfile.electrodes + + for es_name in ("ElectricalSeries1", "ElectricalSeries2"): + es = nwbfile.acquisition[es_name] + assert isinstance(es, ElectricalSeries) + assert es.conversion == 2.34375e-06 + + # make sure the ElectricalSeries objects don't share electrodes + assert not set(nwbfile.acquisition["ElectricalSeries1"].electrodes.data) & set( + nwbfile.acquisition["ElectricalSeries2"].electrodes.data + ) + + assert len(nwbfile.units) == 499 + for col in ("cluster_quality_label", "spike_depths"): + assert col in nwbfile.units + + for es_name in ("ElectricalSeries1", "ElectricalSeries2"): + es = nwbfile.processing["ecephys"].data_interfaces["LFP"][es_name] + assert isinstance(es, ElectricalSeries) + assert es.conversion == 4.6875e-06 + assert es.rate == 2500.0 From 9374d94d1372f6f964ecb3bd628c5f97c18c261b Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 8 Feb 2022 15:40:42 -0400 Subject: [PATCH 29/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 1 - 1 file changed, 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 372ad02c..240f1ea8 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -16,7 +16,6 @@ from tqdm import tqdm from uuid import uuid4 -#from ... import ephys from workflow_array_ephys.pipeline import ephys, probe From 6338dafd0a00c2eb862726977cf6c2470a8c7b1a Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 8 Feb 2022 21:41:00 -0400 Subject: [PATCH 30/77] update import path --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 372ad02c..a62b394d 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -7,7 +7,7 @@ import pynwb from element_session import session from hdmf.backends.hdf5 import H5DataIO -from nwb_conversion_tools.utils.conversion_tools import get_module +from nwb_conversion_tools.utils.nwbfile_tools import get_module from nwb_conversion_tools.utils.genericdatachunkiterator import GenericDataChunkIterator from nwb_conversion_tools.utils.spikeinterfacerecordingdatachunkiterator import ( SpikeInterfaceRecordingDataChunkIterator, From a557b17caaa0d779f91bee105a17fb1418121e00 Mon Sep 17 00:00:00 2001 From: Kabilar Gunalan Date: Mon, 27 Sep 2021 18:20:34 -0500 Subject: [PATCH 31/77] Rebase, squashed. See Details Add element_data_loader for multiple root dirs Update author Fix import Fix OpenEphys session path Update directory path Add print statement Fix for missing `fileTimeSecs` Update error message Suggested adds re upstream components Directing to workflow for upstream `SkullReference` and utility functions --- element_array_ephys/ephys.py | 32 +++++++++++++++++++------------- 1 file changed, 19 insertions(+), 13 deletions(-) diff --git a/element_array_ephys/ephys.py b/element_array_ephys/ephys.py index a434d158..d6efe7d0 100644 --- a/element_array_ephys/ephys.py +++ b/element_array_ephys/ephys.py @@ -57,9 +57,9 @@ def activate(ephys_schema_name, probe_schema_name=None, *, create_schema=True, def get_ephys_root_data_dir() -> list: """ - All data paths, directories in DataJoint Elements are recommended to be - stored as relative paths, with respect to some user-configured "root" - directory, which varies from machine to machine (e.g. different mounted + All data paths, directories in DataJoint Elements are recommended to be + stored as relative paths, with respect to some user-configured "root" + directory, which varies from machine to machine (e.g. different mounted drive locations) get_ephys_root_data_dir() -> list @@ -142,7 +142,7 @@ class EphysFile(dj.Part): def make(self, key): - session_dir = find_full_path(get_ephys_root_data_dir(), + session_dir = find_full_path(get_ephys_root_data_dir(), get_session_directory(key)) inserted_probe_serial_number = (ProbeInsertion * probe.Probe & key).fetch1('probe') @@ -191,7 +191,7 @@ def make(self, key): 'acq_software': acq_software, 'sampling_rate': spikeglx_meta.meta['imSampRate']}) - root_dir = find_root_directory(get_ephys_root_data_dir(), + root_dir = find_root_directory(get_ephys_root_data_dir(), meta_filepath) self.EphysFile.insert1({ **key, @@ -294,8 +294,8 @@ def make(self, key): shank, shank_col, shank_row, _ = spikeglx_recording.apmeta.shankmap['data'][recorded_site] electrode_keys.append(probe_electrodes[(shank, shank_col, shank_row)]) elif acq_software == 'Open Ephys': - - session_dir = find_full_path(get_ephys_root_data_dir(), + + session_dir = find_full_path(get_ephys_root_data_dir(), get_session_directory(key)) loaded_oe = openephys.OpenEphys(session_dir) @@ -457,8 +457,14 @@ class Curation(dj.Manual): def create1_from_clustering_task(self, key, curation_note=''): """ - A function to create a new corresponding "Curation" for a particular +<<<<<<< HEAD + A function to create a new corresponding "Curation" for a particular "ClusteringTask" +======= + A function to create a new corresponding "Curation" for a particular + "ClusteringTask", which assumes that no curation was performed on the + dataset +>>>>>>> dee1c29 (Rebase, squashed. See Details) """ if key not in Clustering(): raise ValueError(f'No corresponding entry in Clustering available' @@ -472,9 +478,9 @@ def create1_from_clustering_task(self, key, curation_note=''): # Synthesize curation_id curation_id = dj.U().aggr(self & key, n='ifnull(max(curation_id)+1,1)').fetch1('n') self.insert1({**key, 'curation_id': curation_id, - 'curation_time': creation_time, + 'curation_time': creation_time, 'curation_output_dir': output_dir, - 'quality_control': is_qc, + 'quality_control': is_qc, 'manual_curation': is_curated, 'curation_note': curation_note}) @@ -622,7 +628,7 @@ def yield_unit_waveforms(): spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) neuropixels_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) elif acq_software == 'Open Ephys': - session_dir = find_full_path(get_ephys_root_data_dir(), + session_dir = find_full_path(get_ephys_root_data_dir(), get_session_directory(key)) openephys_dataset = openephys.OpenEphys(session_dir) neuropixels_recording = openephys_dataset.probes[probe_serial_number] @@ -669,7 +675,7 @@ def get_spikeglx_meta_filepath(ephys_recording_key): except FileNotFoundError: # if not found, search in session_dir again if not spikeglx_meta_filepath.exists(): - session_dir = find_full_path(get_ephys_root_data_dir(), + session_dir = find_full_path(get_ephys_root_data_dir(), get_session_directory( ephys_recording_key)) inserted_probe_serial_number = (ProbeInsertion * probe.Probe @@ -708,7 +714,7 @@ def get_neuropixels_channel2electrode_map(ephys_recording_key, acq_software): for recorded_site, (shank, shank_col, shank_row, _) in enumerate( spikeglx_meta.shankmap['data'])} elif acq_software == 'Open Ephys': - session_dir = find_full_path(get_ephys_root_data_dir(), + session_dir = find_full_path(get_ephys_root_data_dir(), get_session_directory(ephys_recording_key)) openephys_dataset = openephys.OpenEphys(session_dir) probe_serial_number = (ProbeInsertion & ephys_recording_key).fetch1('probe') From 2595daee613f1b080b6bc4a6743865e8f9b42dc7 Mon Sep 17 00:00:00 2001 From: Chris Broz Date: Fri, 11 Feb 2022 16:11:42 -0600 Subject: [PATCH 32/77] Avoid linking_module issues. See details. - add __init__ + schema from git/ttngu207/element_array_ephys@no-curation - add ephys_no_curation schema to match the schema ben has been pulling from - how should we address this in the element? given not currently default - remove unused imports - import datajoint and element_interface.utils find_full_path - add arguments to main export function: - schema names as datajoint config database prefix default - ephys_root_data_dir - default to dj.config or none - add create_virtual_module statements to avoid activate(schema,linking_module=unknown) - declare ephys and probe as global - add assert errors for ephys_root_data_dir!=None when needed - pass ephys_root_data_dir to relevant functions - above permits: from element_array_ephys.export.nwb.nwb import ecephys_session_to_nwb --- element_array_ephys/__init__.py | 12 + element_array_ephys/ephys_no_curation.py | 1024 ++++++++++++++++++++ element_array_ephys/export/nwb/__init__.py | 1 + element_array_ephys/export/nwb/nwb.py | 117 ++- 4 files changed, 1103 insertions(+), 51 deletions(-) create mode 100644 element_array_ephys/ephys_no_curation.py diff --git a/element_array_ephys/__init__.py b/element_array_ephys/__init__.py index e69de29b..c92503b6 100644 --- a/element_array_ephys/__init__.py +++ b/element_array_ephys/__init__.py @@ -0,0 +1,12 @@ +import datajoint as dj +import logging +import os + + +dj.config['enable_python_native_blobs'] = True + + +def get_logger(name): + log = logging.getLogger(name) + log.setLevel(os.getenv('LOGLEVEL', 'INFO')) + return log diff --git a/element_array_ephys/ephys_no_curation.py b/element_array_ephys/ephys_no_curation.py new file mode 100644 index 00000000..0b11645e --- /dev/null +++ b/element_array_ephys/ephys_no_curation.py @@ -0,0 +1,1024 @@ +# pulled from +# https://github.com/ttngu207/element-array-ephys/blob/03cab02ee709e94b151621d00b12e455953dccfb/element_array_ephys/ephys.py + +import datajoint as dj +import pathlib +import re +import numpy as np +import inspect +import importlib +from decimal import Decimal + +from element_interface.utils import find_root_directory, find_full_path, dict_to_uuid + +from .readers import spikeglx, kilosort, openephys +from . import probe, get_logger + + +log = get_logger(__name__) + +schema = dj.schema() + +_linking_module = None + + +def activate(ephys_schema_name, probe_schema_name=None, *, create_schema=True, + create_tables=True, linking_module=None): + """ + activate(ephys_schema_name, probe_schema_name=None, *, create_schema=True, create_tables=True, linking_module=None) + :param ephys_schema_name: schema name on the database server to activate the `ephys` element + :param probe_schema_name: schema name on the database server to activate the `probe` element + - may be omitted if the `probe` element is already activated + :param create_schema: when True (default), create schema in the database if it does not yet exist. + :param create_tables: when True (default), create tables in the database if they do not yet exist. + :param linking_module: a module name or a module containing the + required dependencies to activate the `ephys` element: + Upstream tables: + + Session: table referenced by EphysRecording, typically identifying a recording session + + SkullReference: Reference table for InsertionLocation, specifying the skull reference + used for probe insertion location (e.g. Bregma, Lambda) + Functions: + + get_ephys_root_data_dir() -> list + Retrieve the root data directory - e.g. containing the raw ephys recording files for all subject/sessions. + :return: a string for full path to the root data directory + + get_session_directory(session_key: dict) -> str + Retrieve the session directory containing the recorded Neuropixels data for a given Session + :param session_key: a dictionary of one Session `key` + :return: a string for full path to the session directory + + get_processed_root_data_dir() -> str: + Retrieves the root directory for all processed data to be found from or written to + :return: a string for full path to the root directory for processed data + """ + + if isinstance(linking_module, str): + linking_module = importlib.import_module(linking_module) + assert inspect.ismodule(linking_module),\ + "The argument 'dependency' must be a module's name or a module" + + global _linking_module + _linking_module = linking_module + + probe.activate(probe_schema_name, create_schema=create_schema, + create_tables=create_tables) + schema.activate(ephys_schema_name, create_schema=create_schema, + create_tables=create_tables, add_objects=_linking_module.__dict__) + + +# -------------- Functions required by the elements-ephys --------------- + +def get_ephys_root_data_dir() -> list: + """ + All data paths, directories in DataJoint Elements are recommended to be + stored as relative paths, with respect to some user-configured "root" + directory, which varies from machine to machine (e.g. different mounted + drive locations) + + get_ephys_root_data_dir() -> list + This user-provided function retrieves the possible root data directories + containing the ephys data for all subjects/sessions + (e.g. acquired SpikeGLX or Open Ephys raw files, + output files from spike sorting routines, etc.) + :return: a string for full path to the ephys root data directory, + or list of strings for possible root data directories + """ + root_directories = _linking_module.get_ephys_root_data_dir() + if isinstance(root_directories, (str, pathlib.Path)): + root_directories = [root_directories] + + if hasattr(_linking_module, 'get_processed_root_data_dir'): + root_directories.append(_linking_module.get_processed_root_data_dir()) + + return root_directories + + +def get_session_directory(session_key: dict) -> str: + """ + get_session_directory(session_key: dict) -> str + Retrieve the session directory containing the + recorded Neuropixels data for a given Session + :param session_key: a dictionary of one Session `key` + :return: a string for relative or full path to the session directory + """ + return _linking_module.get_session_directory(session_key) + + +def get_processed_root_data_dir() -> str: + """ + get_processed_root_data_dir() -> str: + Retrieves the root directory for all processed data to be found from or written to + :return: a string for full path to the root directory for processed data + """ + + if hasattr(_linking_module, 'get_processed_root_data_dir'): + return _linking_module.get_processed_root_data_dir() + else: + return get_ephys_root_data_dir()[0] + +# ----------------------------- Table declarations ---------------------- + + +@schema +class AcquisitionSoftware(dj.Lookup): + definition = """ # Name of software used for recording of neuropixels probes - SpikeGLX or Open Ephys + acq_software: varchar(24) + """ + contents = zip(['SpikeGLX', 'Open Ephys']) + + +@schema +class ProbeInsertion(dj.Manual): + definition = """ + # Probe insertion implanted into an animal for a given session. + -> Session + insertion_number: tinyint unsigned + --- + -> probe.Probe + """ + + @classmethod + def auto_generate_entries(cls, session_key): + """ + Method to auto-generate ProbeInsertion entries for a particular session + Probe information is inferred from the meta data found in the session data directory + """ + session_dir = find_full_path(get_ephys_root_data_dir(), + get_session_directory(session_key)) + # search session dir and determine acquisition software + for ephys_pattern, ephys_acq_type in zip(['*.ap.meta', '*.oebin'], + ['SpikeGLX', 'Open Ephys']): + ephys_meta_filepaths = list(session_dir.rglob(ephys_pattern)) + if ephys_meta_filepaths: + acq_software = ephys_acq_type + break + else: + raise FileNotFoundError( + f'Ephys recording data not found!' + f' Neither SpikeGLX nor Open Ephys recording files found in: {session_dir}') + + probe_list, probe_insertion_list = [], [] + if acq_software == 'SpikeGLX': + for meta_fp_idx, meta_filepath in enumerate(ephys_meta_filepaths): + spikeglx_meta = spikeglx.SpikeGLXMeta(meta_filepath) + + probe_key = {'probe_type': spikeglx_meta.probe_model, + 'probe': spikeglx_meta.probe_SN} + if (probe_key['probe'] not in [p['probe'] for p in probe_list] + and probe_key not in probe.Probe()): + probe_list.append(probe_key) + + probe_dir = meta_filepath.parent + try: + probe_number = re.search('(imec)?\d{1}$', probe_dir.name).group() + probe_number = int(probe_number.replace('imec', '')) + except AttributeError: + probe_number = meta_fp_idx + + probe_insertion_list.append({**session_key, + 'probe': spikeglx_meta.probe_SN, + 'insertion_number': int(probe_number)}) + elif acq_software == 'Open Ephys': + loaded_oe = openephys.OpenEphys(session_dir) + for probe_idx, oe_probe in enumerate(loaded_oe.probes.values()): + probe_key = {'probe_type': oe_probe.probe_model, 'probe': oe_probe.probe_SN} + if (probe_key['probe'] not in [p['probe'] for p in probe_list] + and probe_key not in probe.Probe()): + probe_list.append(probe_key) + probe_insertion_list.append({**session_key, + 'probe': oe_probe.probe_SN, + 'insertion_number': probe_idx}) + else: + raise NotImplementedError(f'Unknown acquisition software: {acq_software}') + + probe.Probe.insert(probe_list) + cls.insert(probe_insertion_list, skip_duplicates=True) + + +@schema +class InsertionLocation(dj.Manual): + definition = """ + # Brain Location of a given probe insertion. + -> ProbeInsertion + --- + -> SkullReference + ap_location: decimal(6, 2) # (um) anterior-posterior; ref is 0; more anterior is more positive + ml_location: decimal(6, 2) # (um) medial axis; ref is 0 ; more right is more positive + depth: decimal(6, 2) # (um) manipulator depth relative to surface of the brain (0); more ventral is more negative + theta=null: decimal(5, 2) # (deg) - elevation - rotation about the ml-axis [0, 180] - w.r.t the z+ axis + phi=null: decimal(5, 2) # (deg) - azimuth - rotation about the dv-axis [0, 360] - w.r.t the x+ axis + beta=null: decimal(5, 2) # (deg) rotation about the shank of the probe [-180, 180] - clockwise is increasing in degree - 0 is the probe-front facing anterior + """ + + +@schema +class EphysRecording(dj.Imported): + definition = """ + # Ephys recording from a probe insertion for a given session. + -> ProbeInsertion + --- + -> probe.ElectrodeConfig + -> AcquisitionSoftware + sampling_rate: float # (Hz) + recording_datetime: datetime # datetime of the recording from this probe + recording_duration: float # (seconds) duration of the recording from this probe + """ + + class EphysFile(dj.Part): + definition = """ + # Paths of files of a given EphysRecording round. + -> master + file_path: varchar(255) # filepath relative to root data directory + """ + + def make(self, key): + session_dir = find_full_path(get_ephys_root_data_dir(), + get_session_directory(key)) + inserted_probe_serial_number = (ProbeInsertion * probe.Probe & key).fetch1('probe') + + # search session dir and determine acquisition software + for ephys_pattern, ephys_acq_type in zip(['*.ap.meta', '*.oebin'], + ['SpikeGLX', 'Open Ephys']): + ephys_meta_filepaths = list(session_dir.rglob(ephys_pattern)) + if ephys_meta_filepaths: + acq_software = ephys_acq_type + break + else: + raise FileNotFoundError( + f'Ephys recording data not found!' + f' Neither SpikeGLX nor Open Ephys recording files found' + f' in {session_dir}') + + supported_probe_types = probe.ProbeType.fetch('probe_type') + + if acq_software == 'SpikeGLX': + for meta_filepath in ephys_meta_filepaths: + spikeglx_meta = spikeglx.SpikeGLXMeta(meta_filepath) + if str(spikeglx_meta.probe_SN) == inserted_probe_serial_number: + break + else: + raise FileNotFoundError( + 'No SpikeGLX data found for probe insertion: {}'.format(key)) + + if spikeglx_meta.probe_model in supported_probe_types: + probe_type = spikeglx_meta.probe_model + electrode_query = probe.ProbeType.Electrode & {'probe_type': probe_type} + + probe_electrodes = { + (shank, shank_col, shank_row): key + for key, shank, shank_col, shank_row in zip(*electrode_query.fetch( + 'KEY', 'shank', 'shank_col', 'shank_row'))} + + electrode_group_members = [ + probe_electrodes[(shank, shank_col, shank_row)] + for shank, shank_col, shank_row, _ in spikeglx_meta.shankmap['data']] + else: + raise NotImplementedError( + 'Processing for neuropixels probe model' + ' {} not yet implemented'.format(spikeglx_meta.probe_model)) + + self.insert1({ + **key, + **generate_electrode_config(probe_type, electrode_group_members), + 'acq_software': acq_software, + 'sampling_rate': spikeglx_meta.meta['imSampRate'], + 'recording_datetime': spikeglx_meta.recording_time, + 'recording_duration': (spikeglx_meta.recording_duration + or spikeglx.retrieve_recording_duration(meta_filepath))}) + + root_dir = find_root_directory(get_ephys_root_data_dir(), + meta_filepath) + self.EphysFile.insert1({ + **key, + 'file_path': meta_filepath.relative_to(root_dir).as_posix()}) + elif acq_software == 'Open Ephys': + dataset = openephys.OpenEphys(session_dir) + for serial_number, probe_data in dataset.probes.items(): + if str(serial_number) == inserted_probe_serial_number: + break + else: + raise FileNotFoundError( + 'No Open Ephys data found for probe insertion: {}'.format(key)) + + if probe_data.probe_model in supported_probe_types: + probe_type = probe_data.probe_model + electrode_query = probe.ProbeType.Electrode & {'probe_type': probe_type} + + probe_electrodes = {key['electrode']: key + for key in electrode_query.fetch('KEY')} + + electrode_group_members = [ + probe_electrodes[channel_idx] + for channel_idx in probe_data.ap_meta['channels_indices']] + else: + raise NotImplementedError( + 'Processing for neuropixels' + ' probe model {} not yet implemented'.format(probe_data.probe_model)) + + self.insert1({ + **key, + **generate_electrode_config(probe_type, electrode_group_members), + 'acq_software': acq_software, + 'sampling_rate': probe_data.ap_meta['sample_rate'], + 'recording_datetime': probe_data.recording_info['recording_datetimes'][0], + 'recording_duration': np.sum(probe_data.recording_info['recording_durations'])}) + + root_dir = find_root_directory(get_ephys_root_data_dir(), + probe_data.recording_info['recording_files'][0]) + self.EphysFile.insert([{**key, + 'file_path': fp.relative_to(root_dir).as_posix()} + for fp in probe_data.recording_info['recording_files']]) + else: + raise NotImplementedError(f'Processing ephys files from' + f' acquisition software of type {acq_software} is' + f' not yet implemented') + + +@schema +class LFP(dj.Imported): + definition = """ + # Acquired local field potential (LFP) from a given Ephys recording. + -> EphysRecording + --- + lfp_sampling_rate: float # (Hz) + lfp_time_stamps: longblob # (s) timestamps with respect to the start of the recording (recording_timestamp) + lfp_mean: longblob # (uV) mean of LFP across electrodes - shape (time,) + """ + + class Electrode(dj.Part): + definition = """ + -> master + -> probe.ElectrodeConfig.Electrode + --- + lfp: longblob # (uV) recorded lfp at this electrode + """ + + # Only store LFP for every 9th channel, due to high channel density, + # close-by channels exhibit highly similar LFP + _skip_channel_counts = 9 + + def make(self, key): + acq_software = (EphysRecording * ProbeInsertion & key).fetch1('acq_software') + + electrode_keys, lfp = [], [] + + if acq_software == 'SpikeGLX': + spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) + spikeglx_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) + + lfp_channel_ind = spikeglx_recording.lfmeta.recording_channels[ + -1::-self._skip_channel_counts] + + # Extract LFP data at specified channels and convert to uV + lfp = spikeglx_recording.lf_timeseries[:, lfp_channel_ind] # (sample x channel) + lfp = (lfp * spikeglx_recording.get_channel_bit_volts('lf')[lfp_channel_ind]).T # (channel x sample) + + self.insert1(dict(key, + lfp_sampling_rate=spikeglx_recording.lfmeta.meta['imSampRate'], + lfp_time_stamps=(np.arange(lfp.shape[1]) + / spikeglx_recording.lfmeta.meta['imSampRate']), + lfp_mean=lfp.mean(axis=0))) + + electrode_query = (probe.ProbeType.Electrode + * probe.ElectrodeConfig.Electrode + * EphysRecording & key) + probe_electrodes = { + (shank, shank_col, shank_row): key + for key, shank, shank_col, shank_row in zip(*electrode_query.fetch( + 'KEY', 'shank', 'shank_col', 'shank_row'))} + + for recorded_site in lfp_channel_ind: + shank, shank_col, shank_row, _ = spikeglx_recording.apmeta.shankmap['data'][recorded_site] + electrode_keys.append(probe_electrodes[(shank, shank_col, shank_row)]) + elif acq_software == 'Open Ephys': + oe_probe = get_openephys_probe_data(key) + + lfp_channel_ind = np.r_[ + len(oe_probe.lfp_meta['channels_indices'])-1:0:-self._skip_channel_counts] + + lfp = oe_probe.lfp_timeseries[:, lfp_channel_ind] # (sample x channel) + lfp = (lfp * np.array(oe_probe.lfp_meta['channels_gains'])[lfp_channel_ind]).T # (channel x sample) + lfp_timestamps = oe_probe.lfp_timestamps + + self.insert1(dict(key, + lfp_sampling_rate=oe_probe.lfp_meta['sample_rate'], + lfp_time_stamps=lfp_timestamps, + lfp_mean=lfp.mean(axis=0))) + + electrode_query = (probe.ProbeType.Electrode + * probe.ElectrodeConfig.Electrode + * EphysRecording & key) + probe_electrodes = {key['electrode']: key + for key in electrode_query.fetch('KEY')} + + electrode_keys.extend(probe_electrodes[channel_idx] + for channel_idx in oe_probe.lfp_meta['channels_indices']) + else: + raise NotImplementedError(f'LFP extraction from acquisition software' + f' of type {acq_software} is not yet implemented') + + # single insert in loop to mitigate potential memory issue + for electrode_key, lfp_trace in zip(electrode_keys, lfp): + self.Electrode.insert1({**key, **electrode_key, 'lfp': lfp_trace}) + + +# ------------ Clustering -------------- + +@schema +class ClusteringMethod(dj.Lookup): + definition = """ + # Method for clustering + clustering_method: varchar(16) + --- + clustering_method_desc: varchar(1000) + """ + + contents = [('kilosort2.5', 'kilosort2.5 clustering method'), + ('kilosort3', 'kilosort3 clustering method')] + + +@schema +class ClusteringParamSet(dj.Lookup): + definition = """ + # Parameter set to be used in a clustering procedure + paramset_idx: smallint + --- + -> ClusteringMethod + paramset_desc: varchar(128) + param_set_hash: uuid + unique index (param_set_hash) + params: longblob # dictionary of all applicable parameters + """ + + @classmethod + def insert_new_params(cls, clustering_method: str, paramset_desc: str, + params: dict, paramset_idx: int = None): + if paramset_idx is None: + paramset_idx = (dj.U().aggr(cls, n='max(paramset_idx)').fetch1('n') or 0) + 1 + + param_dict = {'clustering_method': clustering_method, + 'paramset_idx': paramset_idx, + 'paramset_desc': paramset_desc, + 'params': params, + 'param_set_hash': dict_to_uuid( + {**params, 'clustering_method': clustering_method}) + } + param_query = cls & {'param_set_hash': param_dict['param_set_hash']} + + if param_query: # If the specified param-set already exists + existing_paramset_idx = param_query.fetch1('paramset_idx') + if existing_paramset_idx == paramset_idx: # If the existing set has the same paramset_idx: job done + return + else: # If not same name: human error, trying to add the same paramset with different name + raise dj.DataJointError( + f'The specified param-set already exists' + f' - with paramset_idx: {existing_paramset_idx}') + else: + if {'paramset_idx': paramset_idx} in cls.proj(): + raise dj.DataJointError( + f'The specified paramset_idx {paramset_idx} already exists,' + f' please pick a different one.') + cls.insert1(param_dict) + + +@schema +class ClusterQualityLabel(dj.Lookup): + definition = """ + # Quality + cluster_quality_label: varchar(100) # cluster quality type - e.g. 'good', 'MUA', 'noise', etc. + --- + cluster_quality_description: varchar(4000) + """ + contents = [ + ('good', 'single unit'), + ('ok', 'probably a single unit, but could be contaminated'), + ('mua', 'multi-unit activity'), + ('noise', 'bad unit') + ] + + +@schema +class ClusteringTask(dj.Manual): + definition = """ + # Manual table for defining a clustering task ready to be run + -> EphysRecording + -> ClusteringParamSet + --- + clustering_output_dir='': varchar(255) # clustering output directory relative to the clustering root data directory + task_mode='load': enum('load', 'trigger') # 'load': load computed analysis results, 'trigger': trigger computation + """ + + @classmethod + def infer_output_dir(cls, key, relative=False, mkdir=False): + """ + Given a 'key' to an entry in this table + Return the expected clustering_output_dir based on the following convention: + processed_dir / session_dir / probe_{insertion_number} / {clustering_method}_{paramset_idx} + e.g.: sub4/sess1/probe_2/kilosort2_0 + """ + processed_dir = pathlib.Path(get_processed_root_data_dir()) + session_dir = find_full_path(get_ephys_root_data_dir(), + get_session_directory(key)) + root_dir = find_root_directory(get_ephys_root_data_dir(), session_dir) + + method = (ClusteringParamSet * ClusteringMethod & key).fetch1( + 'clustering_method').replace(".", "-") + + output_dir = (processed_dir + / session_dir.relative_to(root_dir) + / f'probe_{key["insertion_number"]}' + / f'{method}_{key["paramset_idx"]}') + + if mkdir: + output_dir.mkdir(parents=True, exist_ok=True) + log.info(f'{output_dir} created!') + + return output_dir.relative_to(processed_dir) if relative else output_dir + + @classmethod + def auto_generate_entries(cls, ephys_recording_key): + """ + Method to auto-generate ClusteringTask entries for a particular ephys recording + Output directory is auto-generated based on the convention + defined in `ClusteringTask.infer_output_dir()` + Default parameter set used: paramset_idx = 0 + """ + key = {**ephys_recording_key, 'paramset_idx': 0} + + processed_dir = get_processed_root_data_dir() + output_dir = ClusteringTask.infer_output_dir(key, relative=False, mkdir=True) + + try: + kilosort.Kilosort(output_dir) # check if the directory is a valid Kilosort output + except FileNotFoundError: + task_mode = 'trigger' + else: + task_mode = 'load' + + cls.insert1({ + **key, + 'clustering_output_dir': output_dir.relative_to(processed_dir).as_posix(), + 'task_mode': task_mode}) + + +@schema +class Clustering(dj.Imported): + """ + A processing table to handle each ClusteringTask: + + If `task_mode == "trigger"`: trigger clustering analysis + according to the ClusteringParamSet (e.g. launch a kilosort job) + + If `task_mode == "load"`: verify output + """ + definition = """ + # Clustering Procedure + -> ClusteringTask + --- + clustering_time: datetime # time of generation of this set of clustering results + package_version='': varchar(16) + """ + + def make(self, key): + task_mode, output_dir = (ClusteringTask & key).fetch1( + 'task_mode', 'clustering_output_dir') + + if not output_dir: + output_dir = ClusteringTask.infer_output_dir(key, relative=True, mkdir=True) + # update clustering_output_dir + ClusteringTask.update1({**key, 'clustering_output_dir': output_dir.as_posix()}) + + kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) + + if task_mode == 'load': + kilosort.Kilosort(kilosort_dir) # check if the directory is a valid Kilosort output + elif task_mode == 'trigger': + acq_software, clustering_method, params = (ClusteringTask * EphysRecording + * ClusteringParamSet & key).fetch1( + 'acq_software', 'clustering_method', 'params') + + if 'kilosort' in clustering_method: + from element_array_ephys.readers import kilosort_triggering + + # add additional probe-recording and channels details into `params` + params = {**params, **get_recording_channels_details(key)} + params['fs'] = params['sample_rate'] + + if acq_software == 'SpikeGLX': + spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) + spikeglx_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) + spikeglx_recording.validate_file('ap') + + if clustering_method.startswith('pykilosort'): + kilosort_triggering.run_pykilosort( + continuous_file=spikeglx_recording.root_dir / ( + spikeglx_recording.root_name + '.ap.bin'), + kilosort_output_directory=kilosort_dir, + channel_ind=params.pop('channel_ind'), + x_coords=params.pop('x_coords'), + y_coords=params.pop('y_coords'), + shank_ind=params.pop('shank_ind'), + connected=params.pop('connected'), + sample_rate=params.pop('sample_rate'), + params=params) + else: + run_kilosort = kilosort_triggering.SGLXKilosortPipeline( + npx_input_dir=spikeglx_meta_filepath.parent, + ks_output_dir=kilosort_dir, + params=params, + KS2ver=f'{Decimal(clustering_method.replace("kilosort", "")):.1f}', + run_CatGT=False) + run_kilosort.run_modules() + elif acq_software == 'Open Ephys': + oe_probe = get_openephys_probe_data(key) + + assert len(oe_probe.recording_info['recording_files']) == 1 + + # run kilosort + if clustering_method.startswith('pykilosort'): + kilosort_triggering.run_pykilosort( + continuous_file=pathlib.Path(oe_probe.recording_info['recording_files'][0]) / 'continuous.dat', + kilosort_output_directory=kilosort_dir, + channel_ind=params.pop('channel_ind'), + x_coords=params.pop('x_coords'), + y_coords=params.pop('y_coords'), + shank_ind=params.pop('shank_ind'), + connected=params.pop('connected'), + sample_rate=params.pop('sample_rate'), + params=params) + else: + run_kilosort = kilosort_triggering.OpenEphysKilosortPipeline( + npx_input_dir=oe_probe.recording_info['recording_files'][0], + ks_output_dir=kilosort_dir, + params=params, + KS2ver=f'{Decimal(clustering_method.replace("kilosort", "")):.1f}') + run_kilosort.run_modules() + else: + raise NotImplementedError(f'Automatic triggering of {clustering_method}' + f' clustering analysis is not yet supported') + + else: + raise ValueError(f'Unknown task mode: {task_mode}') + + creation_time, _, _ = kilosort.extract_clustering_info(kilosort_dir) + self.insert1({**key, 'clustering_time': creation_time}) + + +@schema +class Curation(dj.Manual): + definition = """ + # Manual curation procedure + -> Clustering + curation_id: int + --- + curation_time: datetime # time of generation of this set of curated clustering results + curation_output_dir: varchar(255) # output directory of the curated results, relative to root data directory + quality_control: bool # has this clustering result undergone quality control? + manual_curation: bool # has manual curation been performed on this clustering result? + curation_note='': varchar(2000) + """ + + def create1_from_clustering_task(self, key, curation_note=''): + """ + A function to create a new corresponding "Curation" for a particular + "ClusteringTask" + """ + if key not in Clustering(): + raise ValueError(f'No corresponding entry in Clustering available' + f' for: {key}; do `Clustering.populate(key)`') + + task_mode, output_dir = (ClusteringTask & key).fetch1( + 'task_mode', 'clustering_output_dir') + kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) + + creation_time, is_curated, is_qc = kilosort.extract_clustering_info(kilosort_dir) + # Synthesize curation_id + curation_id = dj.U().aggr(self & key, n='ifnull(max(curation_id)+1,1)').fetch1('n') + self.insert1({**key, 'curation_id': curation_id, + 'curation_time': creation_time, + 'curation_output_dir': output_dir, + 'quality_control': is_qc, + 'manual_curation': is_curated, + 'curation_note': curation_note}) + + +@schema +class CuratedClustering(dj.Imported): + definition = """ + # Clustering results of a curation. + -> Curation + """ + + class Unit(dj.Part): + definition = """ + # Properties of a given unit from a round of clustering (and curation) + -> master + unit: int + --- + -> probe.ElectrodeConfig.Electrode # electrode with highest waveform amplitude for this unit + -> ClusterQualityLabel + spike_count: int # how many spikes in this recording for this unit + spike_times: longblob # (s) spike times of this unit, relative to the start of the EphysRecording + spike_sites : longblob # array of electrode associated with each spike + spike_depths : longblob # (um) array of depths associated with each spike, relative to the (0, 0) of the probe + """ + + def make(self, key): + output_dir = (Curation & key).fetch1('curation_output_dir') + kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) + + kilosort_dataset = kilosort.Kilosort(kilosort_dir) + acq_software, sample_rate = (EphysRecording & key).fetch1( + 'acq_software', 'sampling_rate') + + sample_rate = kilosort_dataset.data['params'].get('sample_rate', sample_rate) + + # ---------- Unit ---------- + # -- Remove 0-spike units + withspike_idx = [i for i, u in enumerate(kilosort_dataset.data['cluster_ids']) + if (kilosort_dataset.data['spike_clusters'] == u).any()] + valid_units = kilosort_dataset.data['cluster_ids'][withspike_idx] + valid_unit_labels = kilosort_dataset.data['cluster_groups'][withspike_idx] + # -- Get channel and electrode-site mapping + channel2electrodes = get_neuropixels_channel2electrode_map(key, acq_software) + + # -- Spike-times -- + # spike_times_sec_adj > spike_times_sec > spike_times + spike_time_key = ('spike_times_sec_adj' if 'spike_times_sec_adj' in kilosort_dataset.data + else 'spike_times_sec' if 'spike_times_sec' + in kilosort_dataset.data else 'spike_times') + spike_times = kilosort_dataset.data[spike_time_key] + kilosort_dataset.extract_spike_depths() + + # -- Spike-sites and Spike-depths -- + spike_sites = np.array([channel2electrodes[s]['electrode'] + for s in kilosort_dataset.data['spike_sites']]) + spike_depths = kilosort_dataset.data['spike_depths'] + + # -- Insert unit, label, peak-chn + units = [] + for unit, unit_lbl in zip(valid_units, valid_unit_labels): + if (kilosort_dataset.data['spike_clusters'] == unit).any(): + unit_channel, _ = kilosort_dataset.get_best_channel(unit) + unit_spike_times = (spike_times[kilosort_dataset.data['spike_clusters'] == unit] + / sample_rate) + spike_count = len(unit_spike_times) + + units.append({ + 'unit': unit, + 'cluster_quality_label': unit_lbl, + **channel2electrodes[unit_channel], + 'spike_times': unit_spike_times, + 'spike_count': spike_count, + 'spike_sites': spike_sites[kilosort_dataset.data['spike_clusters'] == unit], + 'spike_depths': spike_depths[kilosort_dataset.data['spike_clusters'] == unit]}) + + self.insert1(key) + self.Unit.insert([{**key, **u} for u in units]) + + +@schema +class WaveformSet(dj.Imported): + definition = """ + # A set of spike waveforms for units out of a given CuratedClustering + -> CuratedClustering + """ + + class PeakWaveform(dj.Part): + definition = """ + # Mean waveform across spikes for a given unit at its representative electrode + -> master + -> CuratedClustering.Unit + --- + peak_electrode_waveform: longblob # (uV) mean waveform for a given unit at its representative electrode + """ + + class Waveform(dj.Part): + definition = """ + # Spike waveforms and their mean across spikes for the given unit + -> master + -> CuratedClustering.Unit + -> probe.ElectrodeConfig.Electrode + --- + waveform_mean: longblob # (uV) mean waveform across spikes of the given unit + waveforms=null: longblob # (uV) (spike x sample) waveforms of a sampling of spikes at the given electrode for the given unit + """ + + def make(self, key): + output_dir = (Curation & key).fetch1('curation_output_dir') + kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) + + kilosort_dataset = kilosort.Kilosort(kilosort_dir) + + acq_software, probe_serial_number = (EphysRecording * ProbeInsertion & key).fetch1( + 'acq_software', 'probe') + + # -- Get channel and electrode-site mapping + recording_key = (EphysRecording & key).fetch1('KEY') + channel2electrodes = get_neuropixels_channel2electrode_map(recording_key, acq_software) + + is_qc = (Curation & key).fetch1('quality_control') + + # Get all units + units = {u['unit']: u for u in (CuratedClustering.Unit & key).fetch( + as_dict=True, order_by='unit')} + + if is_qc: + unit_waveforms = np.load(kilosort_dir / 'mean_waveforms.npy') # unit x channel x sample + + def yield_unit_waveforms(): + for unit_no, unit_waveform in zip(kilosort_dataset.data['cluster_ids'], + unit_waveforms): + unit_peak_waveform = {} + unit_electrode_waveforms = [] + if unit_no in units: + for channel, channel_waveform in zip( + kilosort_dataset.data['channel_map'], + unit_waveform): + unit_electrode_waveforms.append({ + **units[unit_no], **channel2electrodes[channel], + 'waveform_mean': channel_waveform}) + if channel2electrodes[channel]['electrode'] == units[unit_no]['electrode']: + unit_peak_waveform = { + **units[unit_no], + 'peak_electrode_waveform': channel_waveform} + yield unit_peak_waveform, unit_electrode_waveforms + else: + if acq_software == 'SpikeGLX': + spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) + neuropixels_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) + elif acq_software == 'Open Ephys': + session_dir = find_full_path(get_ephys_root_data_dir(), + get_session_directory(key)) + openephys_dataset = openephys.OpenEphys(session_dir) + neuropixels_recording = openephys_dataset.probes[probe_serial_number] + + def yield_unit_waveforms(): + for unit_dict in units.values(): + unit_peak_waveform = {} + unit_electrode_waveforms = [] + + spikes = unit_dict['spike_times'] + waveforms = neuropixels_recording.extract_spike_waveforms( + spikes, kilosort_dataset.data['channel_map']) # (sample x channel x spike) + waveforms = waveforms.transpose((1, 2, 0)) # (channel x spike x sample) + for channel, channel_waveform in zip( + kilosort_dataset.data['channel_map'], waveforms): + unit_electrode_waveforms.append({ + **unit_dict, **channel2electrodes[channel], + 'waveform_mean': channel_waveform.mean(axis=0), + 'waveforms': channel_waveform}) + if channel2electrodes[channel]['electrode'] == unit_dict['electrode']: + unit_peak_waveform = { + **unit_dict, + 'peak_electrode_waveform': channel_waveform.mean(axis=0)} + + yield unit_peak_waveform, unit_electrode_waveforms + + # insert waveform on a per-unit basis to mitigate potential memory issue + self.insert1(key) + for unit_peak_waveform, unit_electrode_waveforms in yield_unit_waveforms(): + if unit_peak_waveform: + self.PeakWaveform.insert1(unit_peak_waveform, ignore_extra_fields=True) + if unit_electrode_waveforms: + self.Waveform.insert(unit_electrode_waveforms, ignore_extra_fields=True) + + +# ---------------- HELPER FUNCTIONS ---------------- + +def get_spikeglx_meta_filepath(ephys_recording_key): + # attempt to retrieve from EphysRecording.EphysFile + spikeglx_meta_filepath = (EphysRecording.EphysFile & ephys_recording_key + & 'file_path LIKE "%.ap.meta"').fetch1('file_path') + + try: + spikeglx_meta_filepath = find_full_path(get_ephys_root_data_dir(), + spikeglx_meta_filepath) + except FileNotFoundError: + # if not found, search in session_dir again + if not spikeglx_meta_filepath.exists(): + session_dir = find_full_path(get_ephys_root_data_dir(), + get_session_directory( + ephys_recording_key)) + inserted_probe_serial_number = (ProbeInsertion * probe.Probe + & ephys_recording_key).fetch1('probe') + + spikeglx_meta_filepaths = [fp for fp in session_dir.rglob('*.ap.meta')] + for meta_filepath in spikeglx_meta_filepaths: + spikeglx_meta = spikeglx.SpikeGLXMeta(meta_filepath) + if str(spikeglx_meta.probe_SN) == inserted_probe_serial_number: + spikeglx_meta_filepath = meta_filepath + break + else: + raise FileNotFoundError( + 'No SpikeGLX data found for probe insertion: {}'.format(ephys_recording_key)) + + return spikeglx_meta_filepath + + +def get_openephys_probe_data(ephys_recording_key): + inserted_probe_serial_number = (ProbeInsertion * probe.Probe + & ephys_recording_key).fetch1('probe') + session_dir = find_full_path(get_ephys_root_data_dir(), + get_session_directory(ephys_recording_key)) + loaded_oe = openephys.OpenEphys(session_dir) + return loaded_oe.probes[inserted_probe_serial_number] + + +def get_neuropixels_channel2electrode_map(ephys_recording_key, acq_software): + if acq_software == 'SpikeGLX': + spikeglx_meta_filepath = get_spikeglx_meta_filepath(ephys_recording_key) + spikeglx_meta = spikeglx.SpikeGLXMeta(spikeglx_meta_filepath) + electrode_config_key = (EphysRecording * probe.ElectrodeConfig + & ephys_recording_key).fetch1('KEY') + + electrode_query = (probe.ProbeType.Electrode + * probe.ElectrodeConfig.Electrode & electrode_config_key) + + probe_electrodes = { + (shank, shank_col, shank_row): key + for key, shank, shank_col, shank_row in zip(*electrode_query.fetch( + 'KEY', 'shank', 'shank_col', 'shank_row'))} + + channel2electrode_map = { + recorded_site: probe_electrodes[(shank, shank_col, shank_row)] + for recorded_site, (shank, shank_col, shank_row, _) in enumerate( + spikeglx_meta.shankmap['data'])} + elif acq_software == 'Open Ephys': + probe_dataset = get_openephys_probe_data(ephys_recording_key) + + electrode_query = (probe.ProbeType.Electrode + * probe.ElectrodeConfig.Electrode + * EphysRecording & ephys_recording_key) + + probe_electrodes = {key['electrode']: key + for key in electrode_query.fetch('KEY')} + + channel2electrode_map = { + channel_idx: probe_electrodes[channel_idx] + for channel_idx in probe_dataset.ap_meta['channels_indices']} + + return channel2electrode_map + + +def generate_electrode_config(probe_type: str, electrodes: list): + """ + Generate and insert new ElectrodeConfig + :param probe_type: probe type (e.g. neuropixels 2.0 - SS) + :param electrodes: list of the electrode dict (keys of the probe.ProbeType.Electrode table) + :return: a dict representing a key of the probe.ElectrodeConfig table + """ + # compute hash for the electrode config (hash of dict of all ElectrodeConfig.Electrode) + electrode_config_hash = dict_to_uuid({k['electrode']: k for k in electrodes}) + + electrode_list = sorted([k['electrode'] for k in electrodes]) + electrode_gaps = ([-1] + + np.where(np.diff(electrode_list) > 1)[0].tolist() + + [len(electrode_list) - 1]) + electrode_config_name = '; '.join([ + f'{electrode_list[start + 1]}-{electrode_list[end]}' + for start, end in zip(electrode_gaps[:-1], electrode_gaps[1:])]) + + electrode_config_key = {'electrode_config_hash': electrode_config_hash} + + # ---- make new ElectrodeConfig if needed ---- + if not probe.ElectrodeConfig & electrode_config_key: + probe.ElectrodeConfig.insert1({**electrode_config_key, 'probe_type': probe_type, + 'electrode_config_name': electrode_config_name}) + probe.ElectrodeConfig.Electrode.insert({**electrode_config_key, **electrode} + for electrode in electrodes) + + return electrode_config_key + + +def get_recording_channels_details(ephys_recording_key): + channels_details = {} + + acq_software, sample_rate = (EphysRecording & ephys_recording_key).fetch1('acq_software', + 'sampling_rate') + + probe_type = (ProbeInsertion * probe.Probe & ephys_recording_key).fetch1('probe_type') + channels_details['probe_type'] = {'neuropixels 1.0 - 3A': '3A', + 'neuropixels 1.0 - 3B': 'NP1', + 'neuropixels UHD': 'NP1100', + 'neuropixels 2.0 - SS': 'NP21', + 'neuropixels 2.0 - MS': 'NP24'}[probe_type] + + electrode_config_key = (probe.ElectrodeConfig * EphysRecording & ephys_recording_key).fetch1('KEY') + channels_details['channel_ind'], channels_details['x_coords'], channels_details[ + 'y_coords'], channels_details['shank_ind'] = ( + probe.ElectrodeConfig.Electrode * probe.ProbeType.Electrode + & electrode_config_key).fetch('electrode', 'x_coord', 'y_coord', 'shank') + channels_details['sample_rate'] = sample_rate + channels_details['num_channels'] = len(channels_details['channel_ind']) + + if acq_software == 'SpikeGLX': + spikeglx_meta_filepath = get_spikeglx_meta_filepath(ephys_recording_key) + spikeglx_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) + channels_details['uVPerBit'] = spikeglx_recording.get_channel_bit_volts('ap')[0] + channels_details['connected'] = np.array( + [v for *_, v in spikeglx_recording.apmeta.shankmap['data']]) + elif acq_software == 'Open Ephys': + oe_probe = get_openephys_probe_data(ephys_recording_key) + channels_details['uVPerBit'] = oe_probe.ap_meta['channels_gains'][0] + channels_details['connected'] = np.array([ + int(v == 1) for c, v in oe_probe.channels_connected.items() + if c in channels_details['channel_ind']]) + + return channels_details diff --git a/element_array_ephys/export/nwb/__init__.py b/element_array_ephys/export/nwb/__init__.py index e69de29b..e967c7cc 100644 --- a/element_array_ephys/export/nwb/__init__.py +++ b/element_array_ephys/export/nwb/__init__.py @@ -0,0 +1 @@ +from .nwb import ecephys_session_to_nwb, write_nwb diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index b06d25c1..ce954af6 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -1,23 +1,17 @@ import os -import datetime import decimal -import importlib import json import numpy as np import pynwb -from element_session import session +import datajoint as dj +from element_interface.utils import find_full_path from hdmf.backends.hdf5 import H5DataIO from nwb_conversion_tools.utils.nwbfile_tools import get_module from nwb_conversion_tools.utils.genericdatachunkiterator import GenericDataChunkIterator from nwb_conversion_tools.utils.spikeinterfacerecordingdatachunkiterator import ( - SpikeInterfaceRecordingDataChunkIterator, -) + SpikeInterfaceRecordingDataChunkIterator) from spikeinterface import extractors from tqdm import tqdm -from uuid import uuid4 - - -from workflow_array_ephys.pipeline import ephys, probe class DecimalEncoder(json.JSONEncoder): @@ -91,7 +85,6 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): session_key: dict nwbfile: pynwb.NWBFile """ - electrodes_query = probe.ProbeType.Electrode * probe.ElectrodeConfig.Electrode for additional_attribute in ["shank_col", "shank_row", "shank"]: @@ -164,8 +157,7 @@ def create_units_table( units_query, paramset_record, name="units", - desc="data on spiking units", -): + desc="data on spiking units"): """ ephys.CuratedClustering.Unit::unit -> units.id @@ -234,8 +226,7 @@ def create_units_table( def add_ephys_units_to_nwb( - session_key: dict, nwbfile: pynwb.NWBFile, primary_clustering_paramset_idx: int = 0 -): + session_key: dict, nwbfile: pynwb.NWBFile, primary_clustering_paramset_idx: int = 0): """ Add spiking data to NWBFile. @@ -318,8 +309,8 @@ def gains_helper(gains): def add_ephys_recording_to_nwb( - session_key: dict, nwbfile: pynwb.NWBFile, end_frame: int = None -): + session_key: dict, ephys_root_data_dir, + nwbfile: pynwb.NWBFile, end_frame: int = None): """Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically applies lossless compression to the data, to the final file might be smaller than the original, but there is no data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. @@ -348,7 +339,7 @@ def add_ephys_recording_to_nwb( ephys.EphysRecording.EphysFile & ephys_recording_record ).fetch1("file_path") relative_path = relative_path.replace("\\","/") - file_path = ephys.find_full_path(ephys.get_ephys_root_data_dir(), relative_path) + file_path = find_full_path(ephys_root_data_dir, relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": extractor = extractors.read_spikeglx(os.path.split(file_path)[0], "imec.ap") @@ -436,8 +427,7 @@ def add_ephys_lfp_from_dj_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): def add_ephys_lfp_from_source_to_nwb( - session_key: dict, nwbfile: pynwb.NWBFile, end_frame=None -): + session_key: dict, ephys_root_data_dir, nwbfile: pynwb.NWBFile, end_frame=None): """ Read the LFP data directly from the source file. Currently, only works for SpikeGLX data. @@ -472,7 +462,7 @@ def add_ephys_lfp_from_source_to_nwb( ephys.EphysRecording.EphysFile & ephys_recording_record ).fetch1("file_path") relative_path = relative_path.replace("\\","/") - file_path = ephys.find_full_path(ephys.get_ephys_root_data_dir(), relative_path) + file_path = find_full_path(ephys_root_data_dir, relative_path) if ephys_recording_record["acq_software"] == "SpikeGLX": extractor = extractors.read_spikeglx(os.path.split(file_path)[0], "imec.lf") @@ -521,38 +511,61 @@ def ecephys_session_to_nwb( project_key=None, protocol_key=None, nwbfile_kwargs=None, -): + # workflow_module='workflow_array_ephys.pipeline', + # skull_reference='', + lab_schema_name=(dj.config['custom']['database.prefix']+'lab'), + subject_schema_name=(dj.config['custom']['database.prefix']+'subject'), + session_schema_name=(dj.config['custom']['database.prefix']+'session'), + ephys_schema_name=(dj.config['custom']['database.prefix']+'ephys'), + probe_schema_name=(dj.config['custom']['database.prefix']+'probe'), + ephys_root_data_dir=(dj.config.get('custom', {}).get('ephys_root_data_dir', None))): """ - Main function for converting ephys data to NWB + Main function for converting ephys data to NWB - Parameters - ---------- - session_key: dict - subject_id: str - subject_id, used if it cannot be automatically inferred - raw: bool - Whether to include the raw data from source. SpikeGLX and OpenEphys are supported - spikes: bool - Whether to include CuratedClustering - lfp: - "dj" - read LFP data from ephys.LFP - "source" - read LFP data from source (SpikeGLX supported) - False - do not convert LFP - end_frame: int, optional - Used to create small test conversions where large datasets are truncated. - lab_key, project_key, and protocol_key: dictionaries used to look up optional additional metadata - nwbfile_kwargs: dict, optional - - If element-session is not being used, this argument is required and must be a dictionary containing - 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime), - the minimal data for instantiating an NWBFile object. - - - If element-session is being used, this argument can optionally be used to add over overwrite NWBFile fields. + Parameters + ---------- + session_key: dict + subject_id: str + subject_id, used if it cannot be automatically inferred + raw: bool + Whether to include the raw data from source. SpikeGLX and OpenEphys are supported + spikes: bool + Whether to include CuratedClustering + lfp: + "dj" - read LFP data from ephys.LFP + "source" - read LFP data from source (SpikeGLX supported) + False - do not convert LFP + end_frame: int, optional + Used to create small test conversions where large datasets are truncated. + lab_key, project_key, and protocol_key: dictionaries used to look up optional additional metadata + nwbfile_kwargs: dict, optional + - If element-session is not being used, this argument is required and must be a dictionary containing + 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime), + the minimal data for instantiating an NWBFile object. + + - If element-session is being used, this argument can optionally be used to add over overwrite NWBFile fields. """ - if session.schema.is_activated(): + global ephys, probe + ephys = dj.create_virtual_module('ephys',ephys_schema_name) + probe = dj.create_virtual_module('probe',probe_schema_name) + + try: + import element_session + from element_session import session from element_session.export.nwb import session_to_nwb + if not element_session.session.schema.is_activated(): + session = dj.create_virtual_module('session', session_schema_name) + HAVE_ELEMENT_SESSION = True + except ModuleNotFoundError: + HAVE_ELEMENT_SESSION = False + + if HAVE_ELEMENT_SESSION: nwbfile = session_to_nwb( session_key, + lab_schema_name=lab_schema_name, + subject_schema_name=subject_schema_name, + session_schema_name=session_schema_name, subject_id=subject_id, lab_key=lab_key, project_key=project_key, @@ -562,16 +575,16 @@ def ecephys_session_to_nwb( else: if not ( isinstance(nwbfile_kwargs, dict) - and {"session_description", "identifier", "session_start_time"}.issubset(nwbfile_kwargs) - ): + and {"session_description", "identifier", "session_start_time"}.issubset(nwbfile_kwargs)): raise ValueError( "If element-session is not activated, you must include nwbfile_kwargs as a dictionary." - "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)" - ) + "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)") nwbfile = pynwb.NWBFile(**nwbfile_kwargs) if raw: - add_ephys_recording_to_nwb(session_key, nwbfile, end_frame=end_frame) + assert ephys_root_data_dir, 'Need to supply root directory for retrieving raw files' + add_ephys_recording_to_nwb(session_key, ephys_root_data_dir=ephys_root_data_dir, + nwbfile=nwbfile, end_frame=end_frame) if spikes: add_ephys_units_to_nwb(session_key, nwbfile) @@ -580,7 +593,9 @@ def ecephys_session_to_nwb( add_ephys_lfp_from_dj_to_nwb(session_key, nwbfile) if lfp == "source": - add_ephys_lfp_from_source_to_nwb(session_key, nwbfile, end_frame=end_frame) + assert ephys_root_data_dir, 'Need to supply root directory for retrieving raw files' + add_ephys_lfp_from_source_to_nwb(session_key, ephys_root_data_dir=ephys_root_data_dir, + nwbfile=nwbfile, end_frame=end_frame) return nwbfile From bf5e82aa63d771151ddd7348d5d37832cad8ac9d Mon Sep 17 00:00:00 2001 From: Thinh Nguyen Date: Tue, 15 Feb 2022 15:18:21 -0600 Subject: [PATCH 33/77] nwb function specification in linking_module --- element_array_ephys/export/nwb/nwb.py | 55 +++++++++------------------ 1 file changed, 17 insertions(+), 38 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index ce954af6..278d61e1 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -13,6 +13,10 @@ from spikeinterface import extractors from tqdm import tqdm +from .. import probe, ephys + +assert probe.schema.is_activated() and ephys.schema.is_activated() + class DecimalEncoder(json.JSONEncoder): def default(self, o): @@ -502,7 +506,6 @@ def add_ephys_lfp_from_source_to_nwb( def ecephys_session_to_nwb( session_key, - subject_id=None, raw=True, spikes=True, lfp="source", @@ -510,15 +513,7 @@ def ecephys_session_to_nwb( lab_key=None, project_key=None, protocol_key=None, - nwbfile_kwargs=None, - # workflow_module='workflow_array_ephys.pipeline', - # skull_reference='', - lab_schema_name=(dj.config['custom']['database.prefix']+'lab'), - subject_schema_name=(dj.config['custom']['database.prefix']+'subject'), - session_schema_name=(dj.config['custom']['database.prefix']+'session'), - ephys_schema_name=(dj.config['custom']['database.prefix']+'ephys'), - probe_schema_name=(dj.config['custom']['database.prefix']+'probe'), - ephys_root_data_dir=(dj.config.get('custom', {}).get('ephys_root_data_dir', None))): + nwbfile_kwargs=None): """ Main function for converting ephys data to NWB @@ -546,43 +541,28 @@ def ecephys_session_to_nwb( - If element-session is being used, this argument can optionally be used to add over overwrite NWBFile fields. """ - global ephys, probe - ephys = dj.create_virtual_module('ephys',ephys_schema_name) - probe = dj.create_virtual_module('probe',probe_schema_name) - try: - import element_session - from element_session import session - from element_session.export.nwb import session_to_nwb - if not element_session.session.schema.is_activated(): - session = dj.create_virtual_module('session', session_schema_name) - HAVE_ELEMENT_SESSION = True - except ModuleNotFoundError: - HAVE_ELEMENT_SESSION = False - - if HAVE_ELEMENT_SESSION: + session_to_nwb = getattr(ephys._linking_module, 'session_to_nwb') + except AttributeError: + if not ( + isinstance(nwbfile_kwargs, dict) + and {"session_description", "identifier", "session_start_time"}.issubset(nwbfile_kwargs)): + raise ValueError( + "If session_to_nwb is not provided, you must include nwbfile_kwargs as a dictionary." + "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)") + nwbfile = pynwb.NWBFile(**nwbfile_kwargs) + else: nwbfile = session_to_nwb( session_key, - lab_schema_name=lab_schema_name, - subject_schema_name=subject_schema_name, - session_schema_name=session_schema_name, - subject_id=subject_id, lab_key=lab_key, project_key=project_key, protocol_key=protocol_key, additional_nwbfile_kwargs=nwbfile_kwargs, ) - else: - if not ( - isinstance(nwbfile_kwargs, dict) - and {"session_description", "identifier", "session_start_time"}.issubset(nwbfile_kwargs)): - raise ValueError( - "If element-session is not activated, you must include nwbfile_kwargs as a dictionary." - "Required fields are 'session_description' (str), 'identifier' (str), and 'session_start_time' (datetime)") - nwbfile = pynwb.NWBFile(**nwbfile_kwargs) + + ephys_root_data_dir = ephys.get_ephys_root_data_dir() if raw: - assert ephys_root_data_dir, 'Need to supply root directory for retrieving raw files' add_ephys_recording_to_nwb(session_key, ephys_root_data_dir=ephys_root_data_dir, nwbfile=nwbfile, end_frame=end_frame) @@ -593,7 +573,6 @@ def ecephys_session_to_nwb( add_ephys_lfp_from_dj_to_nwb(session_key, nwbfile) if lfp == "source": - assert ephys_root_data_dir, 'Need to supply root directory for retrieving raw files' add_ephys_lfp_from_source_to_nwb(session_key, ephys_root_data_dir=ephys_root_data_dir, nwbfile=nwbfile, end_frame=end_frame) From b4ffe1d68fe068a5bd233ad8c9adb3bda09cf145 Mon Sep 17 00:00:00 2001 From: Thinh Nguyen Date: Wed, 16 Feb 2022 10:08:55 -0600 Subject: [PATCH 34/77] import the correct ephys module that has been activated --- element_array_ephys/__init__.py | 4 ++++ element_array_ephys/{ephys.py => ephys_acute.py} | 0 element_array_ephys/export/nwb/nwb.py | 10 ++++++++-- 3 files changed, 12 insertions(+), 2 deletions(-) rename element_array_ephys/{ephys.py => ephys_acute.py} (100%) diff --git a/element_array_ephys/__init__.py b/element_array_ephys/__init__.py index c92503b6..48b394f3 100644 --- a/element_array_ephys/__init__.py +++ b/element_array_ephys/__init__.py @@ -10,3 +10,7 @@ def get_logger(name): log = logging.getLogger(name) log.setLevel(os.getenv('LOGLEVEL', 'INFO')) return log + + +# ephys_acute as default +import ephys_acute as ephys \ No newline at end of file diff --git a/element_array_ephys/ephys.py b/element_array_ephys/ephys_acute.py similarity index 100% rename from element_array_ephys/ephys.py rename to element_array_ephys/ephys_acute.py diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 278d61e1..943e47e0 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -13,9 +13,15 @@ from spikeinterface import extractors from tqdm import tqdm -from .. import probe, ephys +from .. import probe, ephys_acute, ephys_chronic, ephys_no_curation -assert probe.schema.is_activated() and ephys.schema.is_activated() +assert probe.schema.is_activated(), 'probe not yet activated' + +for ephys in (ephys_acute, ephys_chronic, ephys_no_curation): + if ephys.schema.is_activated(): + break +else: + raise AssertionError('ephys not yet activated') class DecimalEncoder(json.JSONEncoder): From 782a4a5ee0cb7c6e1e382c2985dce45c418c7c19 Mon Sep 17 00:00:00 2001 From: Chris Broz Date: Wed, 16 Feb 2022 11:05:36 -0600 Subject: [PATCH 35/77] adjust imports in __init__ and nwb.py --- element_array_ephys/__init__.py | 2 +- element_array_ephys/export/nwb/nwb.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/element_array_ephys/__init__.py b/element_array_ephys/__init__.py index 48b394f3..aa789eec 100644 --- a/element_array_ephys/__init__.py +++ b/element_array_ephys/__init__.py @@ -13,4 +13,4 @@ def get_logger(name): # ephys_acute as default -import ephys_acute as ephys \ No newline at end of file +import element_array_ephys.ephys_acute as ephys diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 943e47e0..afaf512e 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -13,7 +13,7 @@ from spikeinterface import extractors from tqdm import tqdm -from .. import probe, ephys_acute, ephys_chronic, ephys_no_curation +from ... import probe, ephys_acute, ephys_chronic, ephys_no_curation assert probe.schema.is_activated(), 'probe not yet activated' From fafdde1257aafd00b3c47739be8a73a4c7f09087 Mon Sep 17 00:00:00 2001 From: bendichter Date: Mon, 21 Feb 2022 20:12:27 -0600 Subject: [PATCH 36/77] add tests for getting lfp data from datajoint --- element_array_ephys/export/nwb/tests.py | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/element_array_ephys/export/nwb/tests.py b/element_array_ephys/export/nwb/tests.py index 4525038a..55c088ca 100644 --- a/element_array_ephys/export/nwb/tests.py +++ b/element_array_ephys/export/nwb/tests.py @@ -1,5 +1,6 @@ from element_array_ephys.export.nwb.nwb import ecephys_session_to_nwb, write_nwb +import numpy as np from pynwb.ecephys import ElectricalSeries @@ -35,3 +36,18 @@ def test_convert_to_nwb(): assert isinstance(es, ElectricalSeries) assert es.conversion == 4.6875e-06 assert es.rate == 2500.0 + + +def test_convert_to_nwb_with_dj_lfp(): + nwbfile = ecephys_session_to_nwb( + dict(subject="subject5", session_datetime="2020-05-12 04:13:07"), + lfp="dj", + spikes=False, + ) + + for es_name in ("ElectricalSeries1", "ElectricalSeries2"): + es = nwbfile.processing["ecephys"].data_interfaces["LFP"][es_name] + assert isinstance(es, ElectricalSeries) + assert es.conversion == 1.0 + assert isinstance(es.timestamps, np.ndarray) + From a7f46242d849a65149cb7207a87ec489402f7452 Mon Sep 17 00:00:00 2001 From: bendichter Date: Mon, 28 Feb 2022 15:34:50 -0400 Subject: [PATCH 37/77] import GenericDataChunkIterator from hdmf --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index b06d25c1..b17d4a15 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -7,8 +7,8 @@ import pynwb from element_session import session from hdmf.backends.hdf5 import H5DataIO +from hdmf.data_utils import GenericDataChunkIterator from nwb_conversion_tools.utils.nwbfile_tools import get_module -from nwb_conversion_tools.utils.genericdatachunkiterator import GenericDataChunkIterator from nwb_conversion_tools.utils.spikeinterfacerecordingdatachunkiterator import ( SpikeInterfaceRecordingDataChunkIterator, ) From 806684f4534c11982e78d50264f47257c3c3018f Mon Sep 17 00:00:00 2001 From: bendichter Date: Mon, 28 Feb 2022 15:40:11 -0400 Subject: [PATCH 38/77] rmv subject_id --- element_array_ephys/export/nwb/nwb.py | 4 ---- 1 file changed, 4 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index b17d4a15..22b96209 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -512,7 +512,6 @@ def add_ephys_lfp_from_source_to_nwb( def ecephys_session_to_nwb( session_key, - subject_id=None, raw=True, spikes=True, lfp="source", @@ -528,8 +527,6 @@ def ecephys_session_to_nwb( Parameters ---------- session_key: dict - subject_id: str - subject_id, used if it cannot be automatically inferred raw: bool Whether to include the raw data from source. SpikeGLX and OpenEphys are supported spikes: bool @@ -553,7 +550,6 @@ def ecephys_session_to_nwb( from element_session.export.nwb import session_to_nwb nwbfile = session_to_nwb( session_key, - subject_id=subject_id, lab_key=lab_key, project_key=project_key, protocol_key=protocol_key, From c63c8c42fc640319d4cfc4acf90e1b1595e63e3b Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 1 Mar 2022 13:11:13 -0400 Subject: [PATCH 39/77] rmv tests (they are moved to the ephys workflow) --- element_array_ephys/export/nwb/tests.py | 53 ------------------------- 1 file changed, 53 deletions(-) delete mode 100644 element_array_ephys/export/nwb/tests.py diff --git a/element_array_ephys/export/nwb/tests.py b/element_array_ephys/export/nwb/tests.py deleted file mode 100644 index 55c088ca..00000000 --- a/element_array_ephys/export/nwb/tests.py +++ /dev/null @@ -1,53 +0,0 @@ -from element_array_ephys.export.nwb.nwb import ecephys_session_to_nwb, write_nwb - -import numpy as np -from pynwb.ecephys import ElectricalSeries - - -def test_convert_to_nwb(): - - nwbfile = ecephys_session_to_nwb( - dict(subject="subject5", session_datetime="2020-05-12 04:13:07") - ) - - for x in ("262716621", "714000838"): - assert x in nwbfile.devices - - assert len(nwbfile.electrodes) == 1920 - for col in ("shank", "shank_row", "shank_col"): - assert col in nwbfile.electrodes - - for es_name in ("ElectricalSeries1", "ElectricalSeries2"): - es = nwbfile.acquisition[es_name] - assert isinstance(es, ElectricalSeries) - assert es.conversion == 2.34375e-06 - - # make sure the ElectricalSeries objects don't share electrodes - assert not set(nwbfile.acquisition["ElectricalSeries1"].electrodes.data) & set( - nwbfile.acquisition["ElectricalSeries2"].electrodes.data - ) - - assert len(nwbfile.units) == 499 - for col in ("cluster_quality_label", "spike_depths"): - assert col in nwbfile.units - - for es_name in ("ElectricalSeries1", "ElectricalSeries2"): - es = nwbfile.processing["ecephys"].data_interfaces["LFP"][es_name] - assert isinstance(es, ElectricalSeries) - assert es.conversion == 4.6875e-06 - assert es.rate == 2500.0 - - -def test_convert_to_nwb_with_dj_lfp(): - nwbfile = ecephys_session_to_nwb( - dict(subject="subject5", session_datetime="2020-05-12 04:13:07"), - lfp="dj", - spikes=False, - ) - - for es_name in ("ElectricalSeries1", "ElectricalSeries2"): - es = nwbfile.processing["ecephys"].data_interfaces["LFP"][es_name] - assert isinstance(es, ElectricalSeries) - assert es.conversion == 1.0 - assert isinstance(es.timestamps, np.ndarray) - From cc146716036cf8dba74a241160f0b67ae2c552bd Mon Sep 17 00:00:00 2001 From: bendichter Date: Fri, 4 Mar 2022 10:03:57 -0400 Subject: [PATCH 40/77] fix imports --- element_array_ephys/export/nwb/nwb.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 2f923c60..4e9c796f 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -7,9 +7,10 @@ from element_interface.utils import find_full_path from hdmf.backends.hdf5 import H5DataIO from hdmf.data_utils import GenericDataChunkIterator -from nwb_conversion_tools.utils.nwbfile_tools import get_module -from nwb_conversion_tools.utils.spikeinterfacerecordingdatachunkiterator import ( - SpikeInterfaceRecordingDataChunkIterator) +from nwb_conversion_tools.tools.nwb_helpers import get_module +from nwb_conversion_tools.tools.spikeinterface.spikeinterfacerecordingdatachunkiterator import ( + SpikeInterfaceRecordingDataChunkIterator +) from spikeinterface import extractors from tqdm import tqdm From b373b264080e36e379e775499a1b3166c004b546 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:18:04 -0600 Subject: [PATCH 41/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 4e9c796f..e544e3d1 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -148,7 +148,7 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): group=electrode_group, filtering="unknown", imp=-1.0, - x=np.nan, # to do: populate these values once the CCF element is ready + x=np.nan, y=np.nan, z=np.nan, rel_x=electrode["x_coord"], From e9737435f55f8b3ef48c10ee691c0bfd37dd7e21 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:18:58 -0600 Subject: [PATCH 42/77] Update element_array_ephys/ephys_acute.py Co-authored-by: Chris Brozdowski --- element_array_ephys/ephys_acute.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/element_array_ephys/ephys_acute.py b/element_array_ephys/ephys_acute.py index d6efe7d0..0d614ba6 100644 --- a/element_array_ephys/ephys_acute.py +++ b/element_array_ephys/ephys_acute.py @@ -457,14 +457,9 @@ class Curation(dj.Manual): def create1_from_clustering_task(self, key, curation_note=''): """ -<<<<<<< HEAD - A function to create a new corresponding "Curation" for a particular - "ClusteringTask" -======= A function to create a new corresponding "Curation" for a particular "ClusteringTask", which assumes that no curation was performed on the dataset ->>>>>>> dee1c29 (Rebase, squashed. See Details) """ if key not in Clustering(): raise ValueError(f'No corresponding entry in Clustering available' From 3ae6e2d0f96ac74d77cf5617719abc638c1fba66 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:20:11 -0600 Subject: [PATCH 43/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Chris Brozdowski --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index e544e3d1..705b2078 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -324,7 +324,7 @@ def add_ephys_recording_to_nwb( session_key: dict, ephys_root_data_dir, nwbfile: pynwb.NWBFile, end_frame: int = None): """Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically - applies lossless compression to the data, to the final file might be smaller than the original, but there is no + applies lossless compression to the data, so the final file might be smaller than the original, without data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. source data -> acquisition["ElectricalSeries"] From c68f2cabfe9d362d4a176a80eb8b87f068f8e317 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:20:16 -0600 Subject: [PATCH 44/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 705b2078..a54db7ec 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -395,7 +395,7 @@ def add_ephys_lfp_from_dj_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): """ Read LFP data from the data in element-aray-ephys - ephys.LFP::lfp -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data + ephys.LFP.Electrode::lfp -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data ephys.LFP::lfp_time_stamps -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].timestamps Parameters From 584c738982b53d5988e338302e72f7967b43fe2d Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:20:23 -0600 Subject: [PATCH 45/77] Update element_array_ephys/export/nwb/README.md Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/README.md b/element_array_ephys/export/nwb/README.md index 5b146785..363e9b78 100644 --- a/element_array_ephys/export/nwb/README.md +++ b/element_array_ephys/export/nwb/README.md @@ -44,7 +44,7 @@ which can be called independently: 5. `add_ephys_lfp_from_dj_to_nwb`: Read LFP data from the data in element-array-ephys and convert to NWB. - ephys.LFP::lfp -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data + ephys.LFP.Electrode::lfp -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].data ephys.LFP::lfp_time_stamps -> processing["ecephys"].lfp.electrical_series["ElectricalSeries{insertion_number}"].timestamps 6. `add_ephys_lfp_from_source_to_nwb`: Read the LFP data directly from the source file. Currently, only works for From a1fb1934cb34b1ebedb11a10f514cb1ce24b9e00 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:20:34 -0600 Subject: [PATCH 46/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index a54db7ec..a0f25276 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -587,7 +587,7 @@ def write_nwb(nwbfile, fname, check_read=True): Parameters ---------- nwbfile: pynwb.NWBFile - fname: str + fname: str Absolute path including `*.nwb` extension. check_read: bool If True, PyNWB will try to read the produced NWB file and ensure that it can be read. From 24ac2c038b3e7403a3130e43a385c1bab0acb8f5 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:20:41 -0600 Subject: [PATCH 47/77] Update element_array_ephys/export/nwb/README.md Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/README.md b/element_array_ephys/export/nwb/README.md index 363e9b78..744d1370 100644 --- a/element_array_ephys/export/nwb/README.md +++ b/element_array_ephys/export/nwb/README.md @@ -25,7 +25,7 @@ which can be called independently: probe.ProbeType.Electrode::shank_row -> electrodes["shank_row"] 3. `add_ephys_recording_to_nwb`: Read voltage data directly from source files and iteratively transfer them to the - NWB file. Automatically applies lossless compression to the data, to the final file might be smaller than the original, but there is no + NWB file. Automatically applies lossless compression to the data, so the final file might be smaller than the original, but there is no data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. From abccafa522f930f099eefd49c14e69eb95c56067 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:25:18 -0600 Subject: [PATCH 48/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index a0f25276..7ea07c01 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -325,7 +325,7 @@ def add_ephys_recording_to_nwb( nwbfile: pynwb.NWBFile, end_frame: int = None): """Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically applies lossless compression to the data, so the final file might be smaller than the original, without - data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. + data loss. Currently supports Neuropixels data acquired with SpikeGLX or Open Ephys, and relies on SpikeInterface to read the data. source data -> acquisition["ElectricalSeries"] From 00c369144a50cc34acb0fdab12ea0ac93b6627a4 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Sun, 20 Mar 2022 20:25:25 -0600 Subject: [PATCH 49/77] Update element_array_ephys/export/nwb/README.md Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/README.md b/element_array_ephys/export/nwb/README.md index 744d1370..f6e0068b 100644 --- a/element_array_ephys/export/nwb/README.md +++ b/element_array_ephys/export/nwb/README.md @@ -26,7 +26,7 @@ which can be called independently: 3. `add_ephys_recording_to_nwb`: Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically applies lossless compression to the data, so the final file might be smaller than the original, but there is no - data loss. Currently supports neuropixel and openephys, and relies on SpikeInterface to read the data. + data loss. Currently supports Neuropixels data acquired with SpikeGLX or Open Ephys, and relies on SpikeInterface to read the data. source data -> acquisition["ElectricalSeries"] From 12974ff0df04fb8dc15b650e5cbbd3b51aa6340f Mon Sep 17 00:00:00 2001 From: bendichter Date: Sun, 20 Mar 2022 20:48:03 -0600 Subject: [PATCH 50/77] add docstring for gain_helper --- element_array_ephys/export/nwb/nwb.py | 24 +++++++++++++++++++++++- 1 file changed, 23 insertions(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 7ea07c01..5763d2a6 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -313,6 +313,27 @@ def get_electrodes_mapping(electrodes): def gains_helper(gains): + """ + This handles three different cases for gains: + 1. gains are all 1. In this case, return conversion=1e-6, which applies to all + channels and converts from microvolts to volts. + 2. Gains are all equal, but not 1. In this case, multiply this by 1e-6 to apply this + gain to all channels and convert units to volts. + 3. Gains are different for different channels. In this case use the + `channel_conversion` field in addition to the `conversion` field so that each + channel can be converted to volts using its own individual gain. + + Parameters + ---------- + gains: np.ndarray + + Returns + ------- + dict + conversion : float + channel_conversion : np.ndarray + + """ if all(x == 1 for x in gains): return dict(conversion=1e-6, channel_conversion=None) if all(x == gains[0] for x in gains): @@ -323,7 +344,8 @@ def gains_helper(gains): def add_ephys_recording_to_nwb( session_key: dict, ephys_root_data_dir, nwbfile: pynwb.NWBFile, end_frame: int = None): - """Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically + """ + Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically applies lossless compression to the data, so the final file might be smaller than the original, without data loss. Currently supports Neuropixels data acquired with SpikeGLX or Open Ephys, and relies on SpikeInterface to read the data. From ba6cbcf1c3f03cd6a68b96a5b10e28398d9d60e9 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 10:54:20 -0600 Subject: [PATCH 51/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 5763d2a6..711e0b9b 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -517,7 +517,7 @@ def add_ephys_lfp_from_source_to_nwb( lfp.add_electrical_series( pynwb.ecephys.ElectricalSeries( name=f"ElectricalSeries{ephys_recording_record['insertion_number']}", - description=str(ephys_recording_record), + description=f"LFP from probe {probe_id}", data=SpikeInterfaceRecordingDataChunkIterator(extractor), rate=extractor.get_sampling_frequency(), starting_time=( From 49fac083d5de5c853acb8833460edaa07132638e Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 10:58:04 -0600 Subject: [PATCH 52/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 711e0b9b..0b8f8fa0 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -18,7 +18,7 @@ assert probe.schema.is_activated(), 'probe not yet activated' -for ephys in (ephys_acute, ephys_chronic, ephys_no_curation): +for ephys in (ephys_acute): if ephys.schema.is_activated(): break else: From daccfc4dd7a48142739dfee5b00a3ee7c9624d19 Mon Sep 17 00:00:00 2001 From: bendichter Date: Mon, 21 Mar 2022 11:05:59 -0600 Subject: [PATCH 53/77] specify releases for dependencies --- element_array_ephys/export/nwb/requirements.txt | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/element_array_ephys/export/nwb/requirements.txt b/element_array_ephys/export/nwb/requirements.txt index e4c20c0b..9de1b007 100644 --- a/element_array_ephys/export/nwb/requirements.txt +++ b/element_array_ephys/export/nwb/requirements.txt @@ -1,4 +1,4 @@ dandi -git+https://github.com/catalystneuro/nwb-conversion-tools.git -git+https://github.com/SpikeInterface/spikeinterface.git -git+https://github.com/bendichter/python-neo.git@2d0ba5213f6e2add29dd13b2131650e3325edb23 \ No newline at end of file +nwb-conversion-tools==0.11.1 +spikeinterface==0.93.0 +neo==0.10.2 \ No newline at end of file From 9e72773d701550b9c613dfb208f6543ce227e803 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 16:55:41 -0600 Subject: [PATCH 54/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Chris Brozdowski --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 0b8f8fa0..ad9557a0 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -206,7 +206,7 @@ def create_units_table( ) for unit in tqdm( - (clustering_query @ ephys.CuratedClustering.Unit).fetch(as_dict=True), + (ephys.CuratedClustering.Unit & clustering_query).fetch(as_dict=True), desc=f"creating units table for paramset {paramset_record['paramset_idx']}", ): From 2ee2bd544ed777c537dc18b6421b4764462f3482 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 16:56:34 -0600 Subject: [PATCH 55/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index ad9557a0..adde5193 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -14,7 +14,7 @@ from spikeinterface import extractors from tqdm import tqdm -from ... import probe, ephys_acute, ephys_chronic, ephys_no_curation +from ... import probe, ephys_acute assert probe.schema.is_activated(), 'probe not yet activated' From c200699f35d8f90be07f4a8fb194e33b504d9f78 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 17:10:06 -0600 Subject: [PATCH 56/77] Update element_array_ephys/export/nwb/nwb.py --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index adde5193..ef87aa3a 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -18,7 +18,7 @@ assert probe.schema.is_activated(), 'probe not yet activated' -for ephys in (ephys_acute): +for ephys in (ephys_acute,): if ephys.schema.is_activated(): break else: From acdb5f9d25c4f1c1f6e89dd37fd8e1697327a8e9 Mon Sep 17 00:00:00 2001 From: bendichter Date: Mon, 21 Mar 2022 17:23:09 -0600 Subject: [PATCH 57/77] fix docstring for get_electrodes_mapping --- element_array_ephys/export/nwb/nwb.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 0b8f8fa0..75405bbe 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -293,7 +293,7 @@ def add_ephys_units_to_nwb( def get_electrodes_mapping(electrodes): """ - Create a mapping from the group (shank) and electrode id within that group to the row number of the electrodes + Create a mapping from the probe and electrode id to the row number of the electrodes table. This is used in the construction of the DynamicTableRegion that indicates what rows of the electrodes table correspond to the data in an ElectricalSeries. @@ -342,8 +342,11 @@ def gains_helper(gains): def add_ephys_recording_to_nwb( - session_key: dict, ephys_root_data_dir, - nwbfile: pynwb.NWBFile, end_frame: int = None): + session_key: dict, + ephys_root_data_dir: str, + nwbfile: pynwb.NWBFile, + end_frame: int = None, +): """ Read voltage data directly from source files and iteratively transfer them to the NWB file. Automatically applies lossless compression to the data, so the final file might be smaller than the original, without @@ -354,6 +357,7 @@ def add_ephys_recording_to_nwb( Parameters ---------- session_key: dict + ephys_root_data_dir: str nwbfile: NWBFile end_frame: int, optional Used for small test conversions From cbefcde9f81de6be64eeb31a04631fbb00fe6431 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 17:25:23 -0600 Subject: [PATCH 58/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index dcc09900..d2ed586c 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -113,8 +113,8 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): for this_probe in (ephys.ProbeInsertion * probe.Probe & session_key).fetch( as_dict=True ): - insertion_record = (ephys.InsertionLocation & this_probe).fetch1() - if insertion_record: + insertion_record = (ephys.InsertionLocation & this_probe).fetch() + if len(insertion_record)==1: insert_location = json.dumps( { k: v @@ -123,8 +123,10 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): }, cls=DecimalEncoder, ) - else: + elif len(insertion_record)==0: insert_location = "unknown" + else: + raise DataJointError(f'Found multiple insertion locations for {this_probe}') device = nwbfile.create_device( name=this_probe["probe"], From a01ee9ca8e277265382adb674b707ba67f173c01 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 20:25:12 -0600 Subject: [PATCH 59/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index d2ed586c..0eba5a8e 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -13,7 +13,7 @@ ) from spikeinterface import extractors from tqdm import tqdm - +import warnings from ... import probe, ephys_acute assert probe.schema.is_activated(), 'probe not yet activated' From 2e621c02a4809e1c5a540ddc7e22c07f5fcdec1b Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 21 Mar 2022 20:25:20 -0600 Subject: [PATCH 60/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 0eba5a8e..63a4f83f 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -259,6 +259,7 @@ def add_ephys_units_to_nwb( """ if not ephys.ClusteringTask & session_key: + warnings.warn(f'No unit data exists for session:{session_key}') return if nwbfile.electrodes is None: From 15044cf4c85ae67f9e2a73031980180380a0d974 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 22 Mar 2022 07:57:43 -0600 Subject: [PATCH 61/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 1 - 1 file changed, 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 63a4f83f..f61e5c92 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -265,7 +265,6 @@ def add_ephys_units_to_nwb( if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) - # add additional columns to the units table units_query = ephys.CuratedClustering.Unit() & session_key for paramset_record in ( From 2526812f664e7e58a9db368e293dc8e7927615d7 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 22 Mar 2022 07:58:33 -0600 Subject: [PATCH 62/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index f61e5c92..32048f7a 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -196,6 +196,7 @@ def create_units_table( mapping = get_electrodes_mapping(nwbfile.electrodes) units_table = pynwb.misc.Units(name=name, description=desc) + # add additional columns to the units table for additional_attribute in ["cluster_quality_label", "spike_depths"]: units_table.add_column( name=units_query.heading.attributes[additional_attribute].name, From 766f4eb0962cf36e10290cc1ff4dfb27ad74de87 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 22 Mar 2022 07:59:15 -0600 Subject: [PATCH 63/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 1 - 1 file changed, 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 32048f7a..5fe83624 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -266,7 +266,6 @@ def add_ephys_units_to_nwb( if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) - units_query = ephys.CuratedClustering.Unit() & session_key for paramset_record in ( ephys.ClusteringParamSet & ephys.CuratedClustering & session_key From 842beec80ed5d846aed29a48575db1de1457bdf9 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 22 Mar 2022 07:59:21 -0600 Subject: [PATCH 64/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 1 + 1 file changed, 1 insertion(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 5fe83624..0e2e2493 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -197,6 +197,7 @@ def create_units_table( units_table = pynwb.misc.Units(name=name, description=desc) # add additional columns to the units table + units_query = ephys.CuratedClustering.Unit() & session_key for additional_attribute in ["cluster_quality_label", "spike_depths"]: units_table.add_column( name=units_query.heading.attributes[additional_attribute].name, From 3971fe6632030d789ab6b159936d3f4dc2f5f878 Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 22 Mar 2022 08:28:24 -0600 Subject: [PATCH 65/77] add explanation in docstring of add_ephys_units_to_nwb --- element_array_ephys/export/nwb/nwb.py | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index dcc09900..15bc05fc 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -242,6 +242,17 @@ def add_ephys_units_to_nwb( """ Add spiking data to NWBFile. + In NWB, spiking data is stored in a Units table. The primary Units table is + stored at /units. The spiking data in /units is generally the data used in + downstream analysis. Only a single Units table can be stored at /units. Other Units + tables can be stored in a ProcessingModule at /processing/ecephys/. Any number of + Units tables can be stored in this ProcessingModule as long as they have different + names, and these Units tables can store intermediate processing steps or + alternative curations. + + Use `primary_clustering_paramset_idx` to indicate which clustering is primary. All + others will be stored in /processing/ecephys/. + ephys.CuratedClustering.Unit::unit -> units.id ephys.CuratedClustering.Unit::spike_times -> units["spike_times"] ephys.CuratedClustering.Unit::spike_depths -> units["spike_depths"] From 707adff4e6fbb451e0582115b85b72b36002ba9e Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 22 Mar 2022 08:32:38 -0600 Subject: [PATCH 66/77] add explanation for index parameter --- element_array_ephys/export/nwb/nwb.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index bd7666fc..6e2c24f4 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -199,6 +199,9 @@ def create_units_table( # add additional columns to the units table units_query = ephys.CuratedClustering.Unit() & session_key for additional_attribute in ["cluster_quality_label", "spike_depths"]: + # The `index` parameter indicates whether the column is a "ragged array," i.e. + # whether each row of this column is a vector with potentially different lengths + # for each row. units_table.add_column( name=units_query.heading.attributes[additional_attribute].name, description=units_query.heading.attributes[additional_attribute].comment, From 82c86559e86e64580b740325741a5f62e6cf037f Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 22 Mar 2022 08:49:16 -0600 Subject: [PATCH 67/77] fix insertion record --- element_array_ephys/export/nwb/nwb.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 6e2c24f4..3335e668 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -113,17 +113,17 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): for this_probe in (ephys.ProbeInsertion * probe.Probe & session_key).fetch( as_dict=True ): - insertion_record = (ephys.InsertionLocation & this_probe).fetch() - if len(insertion_record)==1: + insertion_record = (ephys.InsertionLocation & this_probe).fetch(as_dict=True) + if len(insertion_record) == 1: insert_location = json.dumps( { k: v - for k, v in insertion_record.items() + for k, v in insertion_record[0].items() if k not in ephys.InsertionLocation.primary_key }, cls=DecimalEncoder, ) - elif len(insertion_record)==0: + elif len(insertion_record) == 0: insert_location = "unknown" else: raise DataJointError(f'Found multiple insertion locations for {this_probe}') From 0152c5d97dd6b40119312f50da62f31e476dbad4 Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 22 Mar 2022 08:58:50 -0600 Subject: [PATCH 68/77] rmv units_query --- element_array_ephys/export/nwb/nwb.py | 10 ++-------- 1 file changed, 2 insertions(+), 8 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 3335e668..e97c12de 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -167,7 +167,6 @@ def add_electrodes_to_nwb(session_key: dict, nwbfile: pynwb.NWBFile): def create_units_table( session_key: dict, nwbfile: pynwb.NWBFile, - units_query, paramset_record, name="units", desc="data on spiking units"): @@ -184,7 +183,6 @@ def create_units_table( ---------- session_key: dict nwbfile: pynwb.NWBFile - units_query: ephys.CuratedClustering.Unit paramset_record: int name: str, optional default="units" @@ -197,14 +195,13 @@ def create_units_table( units_table = pynwb.misc.Units(name=name, description=desc) # add additional columns to the units table - units_query = ephys.CuratedClustering.Unit() & session_key for additional_attribute in ["cluster_quality_label", "spike_depths"]: # The `index` parameter indicates whether the column is a "ragged array," i.e. # whether each row of this column is a vector with potentially different lengths # for each row. units_table.add_column( - name=units_query.heading.attributes[additional_attribute].name, - description=units_query.heading.attributes[additional_attribute].comment, + name=additional_attribute, + description=ephys.CuratedClustering.Unit.heading.attributes[additional_attribute].comment, index=additional_attribute == "spike_depths", ) @@ -281,7 +278,6 @@ def add_ephys_units_to_nwb( if nwbfile.electrodes is None: add_electrodes_to_nwb(session_key, nwbfile) - for paramset_record in ( ephys.ClusteringParamSet & ephys.CuratedClustering & session_key ).fetch("paramset_idx", "clustering_method", "paramset_desc", as_dict=True): @@ -289,7 +285,6 @@ def add_ephys_units_to_nwb( units_table = create_units_table( session_key, nwbfile, - units_query, paramset_record, desc=paramset_record["paramset_desc"], ) @@ -299,7 +294,6 @@ def add_ephys_units_to_nwb( units_table = create_units_table( session_key, nwbfile, - units_query, paramset_record, name=name, desc=paramset_record["paramset_desc"], From c2004eabf06822d563e06515e50476ac76ec3610 Mon Sep 17 00:00:00 2001 From: bendichter Date: Tue, 22 Mar 2022 09:32:48 -0600 Subject: [PATCH 69/77] trying clustering_query..proj() --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index e97c12de..ea2b7b80 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -210,7 +210,7 @@ def create_units_table( ) for unit in tqdm( - (ephys.CuratedClustering.Unit & clustering_query).fetch(as_dict=True), + (ephys.CuratedClustering.Unit & clustering_query.proj()).fetch(as_dict=True), desc=f"creating units table for paramset {paramset_record['paramset_idx']}", ): From 068ea3d3682cc00d988807e804744102f8c0e359 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 22 Mar 2022 12:06:27 -0600 Subject: [PATCH 70/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index ea2b7b80..160a681c 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -42,7 +42,7 @@ def __init__(self, lfp_electrodes_query, chunk_length: int = 10000): """ Parameters ---------- - lfp_electrodes_query: element_array_ephys.ephys_no_curation.LFP + lfp_electrodes_query: element_array_ephys.ephys.LFP.Electrode chunk_length: int, optional Chunks are blocks of disk space where data are stored contiguously and compressed """ From 59028cb22c1ee533d4635d677eea897a933cbf71 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Tue, 22 Mar 2022 13:29:24 -0600 Subject: [PATCH 71/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 160a681c..a858ccd7 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -51,7 +51,7 @@ def __init__(self, lfp_electrodes_query, chunk_length: int = 10000): first_record = ( self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) - ).fetch(as_dict=True)[0] + ).fetch1(as_dict=True) self.n_channels = len(lfp_electrodes_query) self.n_tt = len(first_record["lfp"]) From d94453b1e7704ef41f830ae4c7e06f569a37f545 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 18 Apr 2022 11:52:32 -0400 Subject: [PATCH 72/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index a858ccd7..332fd4ee 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -215,10 +215,10 @@ def create_units_table( ): probe_id, shank_num = ( - probe.ProbeType.Electrode - * ephys.ProbeInsertion + ephys.ProbeInsertion * ephys.CuratedClustering.Unit - & {"unit": unit["unit"], "insertion_number": unit["insertion_number"]} + * probe.ProbeType.Electrode + & unit ).fetch1("probe", "shank") waveform_mean = ( From fb983274e8294fbee5703e6ccaaa3d46ad1394b4 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 18 Apr 2022 11:52:46 -0400 Subject: [PATCH 73/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 332fd4ee..b47a88d1 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -14,7 +14,7 @@ from spikeinterface import extractors from tqdm import tqdm import warnings -from ... import probe, ephys_acute +from ... import probe, ephys_no_curation assert probe.schema.is_activated(), 'probe not yet activated' From 9f9872c37eb325b441703d4204889f6298d1ba4e Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Mon, 18 Apr 2022 11:52:56 -0400 Subject: [PATCH 74/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Kabilar Gunalan --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index b47a88d1..c0a8c8a7 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -18,7 +18,7 @@ assert probe.schema.is_activated(), 'probe not yet activated' -for ephys in (ephys_acute,): +for ephys in (ephys_no_curation,): if ephys.schema.is_activated(): break else: From 3e07c61b7556fffea3c7daa7409f70f51541e76e Mon Sep 17 00:00:00 2001 From: bendichter Date: Wed, 27 Apr 2022 10:44:55 -0400 Subject: [PATCH 75/77] remove ephys_no_curation.py --- element_array_ephys/ephys_no_curation.py | 1024 ---------------------- 1 file changed, 1024 deletions(-) delete mode 100644 element_array_ephys/ephys_no_curation.py diff --git a/element_array_ephys/ephys_no_curation.py b/element_array_ephys/ephys_no_curation.py deleted file mode 100644 index 0b11645e..00000000 --- a/element_array_ephys/ephys_no_curation.py +++ /dev/null @@ -1,1024 +0,0 @@ -# pulled from -# https://github.com/ttngu207/element-array-ephys/blob/03cab02ee709e94b151621d00b12e455953dccfb/element_array_ephys/ephys.py - -import datajoint as dj -import pathlib -import re -import numpy as np -import inspect -import importlib -from decimal import Decimal - -from element_interface.utils import find_root_directory, find_full_path, dict_to_uuid - -from .readers import spikeglx, kilosort, openephys -from . import probe, get_logger - - -log = get_logger(__name__) - -schema = dj.schema() - -_linking_module = None - - -def activate(ephys_schema_name, probe_schema_name=None, *, create_schema=True, - create_tables=True, linking_module=None): - """ - activate(ephys_schema_name, probe_schema_name=None, *, create_schema=True, create_tables=True, linking_module=None) - :param ephys_schema_name: schema name on the database server to activate the `ephys` element - :param probe_schema_name: schema name on the database server to activate the `probe` element - - may be omitted if the `probe` element is already activated - :param create_schema: when True (default), create schema in the database if it does not yet exist. - :param create_tables: when True (default), create tables in the database if they do not yet exist. - :param linking_module: a module name or a module containing the - required dependencies to activate the `ephys` element: - Upstream tables: - + Session: table referenced by EphysRecording, typically identifying a recording session - + SkullReference: Reference table for InsertionLocation, specifying the skull reference - used for probe insertion location (e.g. Bregma, Lambda) - Functions: - + get_ephys_root_data_dir() -> list - Retrieve the root data directory - e.g. containing the raw ephys recording files for all subject/sessions. - :return: a string for full path to the root data directory - + get_session_directory(session_key: dict) -> str - Retrieve the session directory containing the recorded Neuropixels data for a given Session - :param session_key: a dictionary of one Session `key` - :return: a string for full path to the session directory - + get_processed_root_data_dir() -> str: - Retrieves the root directory for all processed data to be found from or written to - :return: a string for full path to the root directory for processed data - """ - - if isinstance(linking_module, str): - linking_module = importlib.import_module(linking_module) - assert inspect.ismodule(linking_module),\ - "The argument 'dependency' must be a module's name or a module" - - global _linking_module - _linking_module = linking_module - - probe.activate(probe_schema_name, create_schema=create_schema, - create_tables=create_tables) - schema.activate(ephys_schema_name, create_schema=create_schema, - create_tables=create_tables, add_objects=_linking_module.__dict__) - - -# -------------- Functions required by the elements-ephys --------------- - -def get_ephys_root_data_dir() -> list: - """ - All data paths, directories in DataJoint Elements are recommended to be - stored as relative paths, with respect to some user-configured "root" - directory, which varies from machine to machine (e.g. different mounted - drive locations) - - get_ephys_root_data_dir() -> list - This user-provided function retrieves the possible root data directories - containing the ephys data for all subjects/sessions - (e.g. acquired SpikeGLX or Open Ephys raw files, - output files from spike sorting routines, etc.) - :return: a string for full path to the ephys root data directory, - or list of strings for possible root data directories - """ - root_directories = _linking_module.get_ephys_root_data_dir() - if isinstance(root_directories, (str, pathlib.Path)): - root_directories = [root_directories] - - if hasattr(_linking_module, 'get_processed_root_data_dir'): - root_directories.append(_linking_module.get_processed_root_data_dir()) - - return root_directories - - -def get_session_directory(session_key: dict) -> str: - """ - get_session_directory(session_key: dict) -> str - Retrieve the session directory containing the - recorded Neuropixels data for a given Session - :param session_key: a dictionary of one Session `key` - :return: a string for relative or full path to the session directory - """ - return _linking_module.get_session_directory(session_key) - - -def get_processed_root_data_dir() -> str: - """ - get_processed_root_data_dir() -> str: - Retrieves the root directory for all processed data to be found from or written to - :return: a string for full path to the root directory for processed data - """ - - if hasattr(_linking_module, 'get_processed_root_data_dir'): - return _linking_module.get_processed_root_data_dir() - else: - return get_ephys_root_data_dir()[0] - -# ----------------------------- Table declarations ---------------------- - - -@schema -class AcquisitionSoftware(dj.Lookup): - definition = """ # Name of software used for recording of neuropixels probes - SpikeGLX or Open Ephys - acq_software: varchar(24) - """ - contents = zip(['SpikeGLX', 'Open Ephys']) - - -@schema -class ProbeInsertion(dj.Manual): - definition = """ - # Probe insertion implanted into an animal for a given session. - -> Session - insertion_number: tinyint unsigned - --- - -> probe.Probe - """ - - @classmethod - def auto_generate_entries(cls, session_key): - """ - Method to auto-generate ProbeInsertion entries for a particular session - Probe information is inferred from the meta data found in the session data directory - """ - session_dir = find_full_path(get_ephys_root_data_dir(), - get_session_directory(session_key)) - # search session dir and determine acquisition software - for ephys_pattern, ephys_acq_type in zip(['*.ap.meta', '*.oebin'], - ['SpikeGLX', 'Open Ephys']): - ephys_meta_filepaths = list(session_dir.rglob(ephys_pattern)) - if ephys_meta_filepaths: - acq_software = ephys_acq_type - break - else: - raise FileNotFoundError( - f'Ephys recording data not found!' - f' Neither SpikeGLX nor Open Ephys recording files found in: {session_dir}') - - probe_list, probe_insertion_list = [], [] - if acq_software == 'SpikeGLX': - for meta_fp_idx, meta_filepath in enumerate(ephys_meta_filepaths): - spikeglx_meta = spikeglx.SpikeGLXMeta(meta_filepath) - - probe_key = {'probe_type': spikeglx_meta.probe_model, - 'probe': spikeglx_meta.probe_SN} - if (probe_key['probe'] not in [p['probe'] for p in probe_list] - and probe_key not in probe.Probe()): - probe_list.append(probe_key) - - probe_dir = meta_filepath.parent - try: - probe_number = re.search('(imec)?\d{1}$', probe_dir.name).group() - probe_number = int(probe_number.replace('imec', '')) - except AttributeError: - probe_number = meta_fp_idx - - probe_insertion_list.append({**session_key, - 'probe': spikeglx_meta.probe_SN, - 'insertion_number': int(probe_number)}) - elif acq_software == 'Open Ephys': - loaded_oe = openephys.OpenEphys(session_dir) - for probe_idx, oe_probe in enumerate(loaded_oe.probes.values()): - probe_key = {'probe_type': oe_probe.probe_model, 'probe': oe_probe.probe_SN} - if (probe_key['probe'] not in [p['probe'] for p in probe_list] - and probe_key not in probe.Probe()): - probe_list.append(probe_key) - probe_insertion_list.append({**session_key, - 'probe': oe_probe.probe_SN, - 'insertion_number': probe_idx}) - else: - raise NotImplementedError(f'Unknown acquisition software: {acq_software}') - - probe.Probe.insert(probe_list) - cls.insert(probe_insertion_list, skip_duplicates=True) - - -@schema -class InsertionLocation(dj.Manual): - definition = """ - # Brain Location of a given probe insertion. - -> ProbeInsertion - --- - -> SkullReference - ap_location: decimal(6, 2) # (um) anterior-posterior; ref is 0; more anterior is more positive - ml_location: decimal(6, 2) # (um) medial axis; ref is 0 ; more right is more positive - depth: decimal(6, 2) # (um) manipulator depth relative to surface of the brain (0); more ventral is more negative - theta=null: decimal(5, 2) # (deg) - elevation - rotation about the ml-axis [0, 180] - w.r.t the z+ axis - phi=null: decimal(5, 2) # (deg) - azimuth - rotation about the dv-axis [0, 360] - w.r.t the x+ axis - beta=null: decimal(5, 2) # (deg) rotation about the shank of the probe [-180, 180] - clockwise is increasing in degree - 0 is the probe-front facing anterior - """ - - -@schema -class EphysRecording(dj.Imported): - definition = """ - # Ephys recording from a probe insertion for a given session. - -> ProbeInsertion - --- - -> probe.ElectrodeConfig - -> AcquisitionSoftware - sampling_rate: float # (Hz) - recording_datetime: datetime # datetime of the recording from this probe - recording_duration: float # (seconds) duration of the recording from this probe - """ - - class EphysFile(dj.Part): - definition = """ - # Paths of files of a given EphysRecording round. - -> master - file_path: varchar(255) # filepath relative to root data directory - """ - - def make(self, key): - session_dir = find_full_path(get_ephys_root_data_dir(), - get_session_directory(key)) - inserted_probe_serial_number = (ProbeInsertion * probe.Probe & key).fetch1('probe') - - # search session dir and determine acquisition software - for ephys_pattern, ephys_acq_type in zip(['*.ap.meta', '*.oebin'], - ['SpikeGLX', 'Open Ephys']): - ephys_meta_filepaths = list(session_dir.rglob(ephys_pattern)) - if ephys_meta_filepaths: - acq_software = ephys_acq_type - break - else: - raise FileNotFoundError( - f'Ephys recording data not found!' - f' Neither SpikeGLX nor Open Ephys recording files found' - f' in {session_dir}') - - supported_probe_types = probe.ProbeType.fetch('probe_type') - - if acq_software == 'SpikeGLX': - for meta_filepath in ephys_meta_filepaths: - spikeglx_meta = spikeglx.SpikeGLXMeta(meta_filepath) - if str(spikeglx_meta.probe_SN) == inserted_probe_serial_number: - break - else: - raise FileNotFoundError( - 'No SpikeGLX data found for probe insertion: {}'.format(key)) - - if spikeglx_meta.probe_model in supported_probe_types: - probe_type = spikeglx_meta.probe_model - electrode_query = probe.ProbeType.Electrode & {'probe_type': probe_type} - - probe_electrodes = { - (shank, shank_col, shank_row): key - for key, shank, shank_col, shank_row in zip(*electrode_query.fetch( - 'KEY', 'shank', 'shank_col', 'shank_row'))} - - electrode_group_members = [ - probe_electrodes[(shank, shank_col, shank_row)] - for shank, shank_col, shank_row, _ in spikeglx_meta.shankmap['data']] - else: - raise NotImplementedError( - 'Processing for neuropixels probe model' - ' {} not yet implemented'.format(spikeglx_meta.probe_model)) - - self.insert1({ - **key, - **generate_electrode_config(probe_type, electrode_group_members), - 'acq_software': acq_software, - 'sampling_rate': spikeglx_meta.meta['imSampRate'], - 'recording_datetime': spikeglx_meta.recording_time, - 'recording_duration': (spikeglx_meta.recording_duration - or spikeglx.retrieve_recording_duration(meta_filepath))}) - - root_dir = find_root_directory(get_ephys_root_data_dir(), - meta_filepath) - self.EphysFile.insert1({ - **key, - 'file_path': meta_filepath.relative_to(root_dir).as_posix()}) - elif acq_software == 'Open Ephys': - dataset = openephys.OpenEphys(session_dir) - for serial_number, probe_data in dataset.probes.items(): - if str(serial_number) == inserted_probe_serial_number: - break - else: - raise FileNotFoundError( - 'No Open Ephys data found for probe insertion: {}'.format(key)) - - if probe_data.probe_model in supported_probe_types: - probe_type = probe_data.probe_model - electrode_query = probe.ProbeType.Electrode & {'probe_type': probe_type} - - probe_electrodes = {key['electrode']: key - for key in electrode_query.fetch('KEY')} - - electrode_group_members = [ - probe_electrodes[channel_idx] - for channel_idx in probe_data.ap_meta['channels_indices']] - else: - raise NotImplementedError( - 'Processing for neuropixels' - ' probe model {} not yet implemented'.format(probe_data.probe_model)) - - self.insert1({ - **key, - **generate_electrode_config(probe_type, electrode_group_members), - 'acq_software': acq_software, - 'sampling_rate': probe_data.ap_meta['sample_rate'], - 'recording_datetime': probe_data.recording_info['recording_datetimes'][0], - 'recording_duration': np.sum(probe_data.recording_info['recording_durations'])}) - - root_dir = find_root_directory(get_ephys_root_data_dir(), - probe_data.recording_info['recording_files'][0]) - self.EphysFile.insert([{**key, - 'file_path': fp.relative_to(root_dir).as_posix()} - for fp in probe_data.recording_info['recording_files']]) - else: - raise NotImplementedError(f'Processing ephys files from' - f' acquisition software of type {acq_software} is' - f' not yet implemented') - - -@schema -class LFP(dj.Imported): - definition = """ - # Acquired local field potential (LFP) from a given Ephys recording. - -> EphysRecording - --- - lfp_sampling_rate: float # (Hz) - lfp_time_stamps: longblob # (s) timestamps with respect to the start of the recording (recording_timestamp) - lfp_mean: longblob # (uV) mean of LFP across electrodes - shape (time,) - """ - - class Electrode(dj.Part): - definition = """ - -> master - -> probe.ElectrodeConfig.Electrode - --- - lfp: longblob # (uV) recorded lfp at this electrode - """ - - # Only store LFP for every 9th channel, due to high channel density, - # close-by channels exhibit highly similar LFP - _skip_channel_counts = 9 - - def make(self, key): - acq_software = (EphysRecording * ProbeInsertion & key).fetch1('acq_software') - - electrode_keys, lfp = [], [] - - if acq_software == 'SpikeGLX': - spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) - spikeglx_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) - - lfp_channel_ind = spikeglx_recording.lfmeta.recording_channels[ - -1::-self._skip_channel_counts] - - # Extract LFP data at specified channels and convert to uV - lfp = spikeglx_recording.lf_timeseries[:, lfp_channel_ind] # (sample x channel) - lfp = (lfp * spikeglx_recording.get_channel_bit_volts('lf')[lfp_channel_ind]).T # (channel x sample) - - self.insert1(dict(key, - lfp_sampling_rate=spikeglx_recording.lfmeta.meta['imSampRate'], - lfp_time_stamps=(np.arange(lfp.shape[1]) - / spikeglx_recording.lfmeta.meta['imSampRate']), - lfp_mean=lfp.mean(axis=0))) - - electrode_query = (probe.ProbeType.Electrode - * probe.ElectrodeConfig.Electrode - * EphysRecording & key) - probe_electrodes = { - (shank, shank_col, shank_row): key - for key, shank, shank_col, shank_row in zip(*electrode_query.fetch( - 'KEY', 'shank', 'shank_col', 'shank_row'))} - - for recorded_site in lfp_channel_ind: - shank, shank_col, shank_row, _ = spikeglx_recording.apmeta.shankmap['data'][recorded_site] - electrode_keys.append(probe_electrodes[(shank, shank_col, shank_row)]) - elif acq_software == 'Open Ephys': - oe_probe = get_openephys_probe_data(key) - - lfp_channel_ind = np.r_[ - len(oe_probe.lfp_meta['channels_indices'])-1:0:-self._skip_channel_counts] - - lfp = oe_probe.lfp_timeseries[:, lfp_channel_ind] # (sample x channel) - lfp = (lfp * np.array(oe_probe.lfp_meta['channels_gains'])[lfp_channel_ind]).T # (channel x sample) - lfp_timestamps = oe_probe.lfp_timestamps - - self.insert1(dict(key, - lfp_sampling_rate=oe_probe.lfp_meta['sample_rate'], - lfp_time_stamps=lfp_timestamps, - lfp_mean=lfp.mean(axis=0))) - - electrode_query = (probe.ProbeType.Electrode - * probe.ElectrodeConfig.Electrode - * EphysRecording & key) - probe_electrodes = {key['electrode']: key - for key in electrode_query.fetch('KEY')} - - electrode_keys.extend(probe_electrodes[channel_idx] - for channel_idx in oe_probe.lfp_meta['channels_indices']) - else: - raise NotImplementedError(f'LFP extraction from acquisition software' - f' of type {acq_software} is not yet implemented') - - # single insert in loop to mitigate potential memory issue - for electrode_key, lfp_trace in zip(electrode_keys, lfp): - self.Electrode.insert1({**key, **electrode_key, 'lfp': lfp_trace}) - - -# ------------ Clustering -------------- - -@schema -class ClusteringMethod(dj.Lookup): - definition = """ - # Method for clustering - clustering_method: varchar(16) - --- - clustering_method_desc: varchar(1000) - """ - - contents = [('kilosort2.5', 'kilosort2.5 clustering method'), - ('kilosort3', 'kilosort3 clustering method')] - - -@schema -class ClusteringParamSet(dj.Lookup): - definition = """ - # Parameter set to be used in a clustering procedure - paramset_idx: smallint - --- - -> ClusteringMethod - paramset_desc: varchar(128) - param_set_hash: uuid - unique index (param_set_hash) - params: longblob # dictionary of all applicable parameters - """ - - @classmethod - def insert_new_params(cls, clustering_method: str, paramset_desc: str, - params: dict, paramset_idx: int = None): - if paramset_idx is None: - paramset_idx = (dj.U().aggr(cls, n='max(paramset_idx)').fetch1('n') or 0) + 1 - - param_dict = {'clustering_method': clustering_method, - 'paramset_idx': paramset_idx, - 'paramset_desc': paramset_desc, - 'params': params, - 'param_set_hash': dict_to_uuid( - {**params, 'clustering_method': clustering_method}) - } - param_query = cls & {'param_set_hash': param_dict['param_set_hash']} - - if param_query: # If the specified param-set already exists - existing_paramset_idx = param_query.fetch1('paramset_idx') - if existing_paramset_idx == paramset_idx: # If the existing set has the same paramset_idx: job done - return - else: # If not same name: human error, trying to add the same paramset with different name - raise dj.DataJointError( - f'The specified param-set already exists' - f' - with paramset_idx: {existing_paramset_idx}') - else: - if {'paramset_idx': paramset_idx} in cls.proj(): - raise dj.DataJointError( - f'The specified paramset_idx {paramset_idx} already exists,' - f' please pick a different one.') - cls.insert1(param_dict) - - -@schema -class ClusterQualityLabel(dj.Lookup): - definition = """ - # Quality - cluster_quality_label: varchar(100) # cluster quality type - e.g. 'good', 'MUA', 'noise', etc. - --- - cluster_quality_description: varchar(4000) - """ - contents = [ - ('good', 'single unit'), - ('ok', 'probably a single unit, but could be contaminated'), - ('mua', 'multi-unit activity'), - ('noise', 'bad unit') - ] - - -@schema -class ClusteringTask(dj.Manual): - definition = """ - # Manual table for defining a clustering task ready to be run - -> EphysRecording - -> ClusteringParamSet - --- - clustering_output_dir='': varchar(255) # clustering output directory relative to the clustering root data directory - task_mode='load': enum('load', 'trigger') # 'load': load computed analysis results, 'trigger': trigger computation - """ - - @classmethod - def infer_output_dir(cls, key, relative=False, mkdir=False): - """ - Given a 'key' to an entry in this table - Return the expected clustering_output_dir based on the following convention: - processed_dir / session_dir / probe_{insertion_number} / {clustering_method}_{paramset_idx} - e.g.: sub4/sess1/probe_2/kilosort2_0 - """ - processed_dir = pathlib.Path(get_processed_root_data_dir()) - session_dir = find_full_path(get_ephys_root_data_dir(), - get_session_directory(key)) - root_dir = find_root_directory(get_ephys_root_data_dir(), session_dir) - - method = (ClusteringParamSet * ClusteringMethod & key).fetch1( - 'clustering_method').replace(".", "-") - - output_dir = (processed_dir - / session_dir.relative_to(root_dir) - / f'probe_{key["insertion_number"]}' - / f'{method}_{key["paramset_idx"]}') - - if mkdir: - output_dir.mkdir(parents=True, exist_ok=True) - log.info(f'{output_dir} created!') - - return output_dir.relative_to(processed_dir) if relative else output_dir - - @classmethod - def auto_generate_entries(cls, ephys_recording_key): - """ - Method to auto-generate ClusteringTask entries for a particular ephys recording - Output directory is auto-generated based on the convention - defined in `ClusteringTask.infer_output_dir()` - Default parameter set used: paramset_idx = 0 - """ - key = {**ephys_recording_key, 'paramset_idx': 0} - - processed_dir = get_processed_root_data_dir() - output_dir = ClusteringTask.infer_output_dir(key, relative=False, mkdir=True) - - try: - kilosort.Kilosort(output_dir) # check if the directory is a valid Kilosort output - except FileNotFoundError: - task_mode = 'trigger' - else: - task_mode = 'load' - - cls.insert1({ - **key, - 'clustering_output_dir': output_dir.relative_to(processed_dir).as_posix(), - 'task_mode': task_mode}) - - -@schema -class Clustering(dj.Imported): - """ - A processing table to handle each ClusteringTask: - + If `task_mode == "trigger"`: trigger clustering analysis - according to the ClusteringParamSet (e.g. launch a kilosort job) - + If `task_mode == "load"`: verify output - """ - definition = """ - # Clustering Procedure - -> ClusteringTask - --- - clustering_time: datetime # time of generation of this set of clustering results - package_version='': varchar(16) - """ - - def make(self, key): - task_mode, output_dir = (ClusteringTask & key).fetch1( - 'task_mode', 'clustering_output_dir') - - if not output_dir: - output_dir = ClusteringTask.infer_output_dir(key, relative=True, mkdir=True) - # update clustering_output_dir - ClusteringTask.update1({**key, 'clustering_output_dir': output_dir.as_posix()}) - - kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) - - if task_mode == 'load': - kilosort.Kilosort(kilosort_dir) # check if the directory is a valid Kilosort output - elif task_mode == 'trigger': - acq_software, clustering_method, params = (ClusteringTask * EphysRecording - * ClusteringParamSet & key).fetch1( - 'acq_software', 'clustering_method', 'params') - - if 'kilosort' in clustering_method: - from element_array_ephys.readers import kilosort_triggering - - # add additional probe-recording and channels details into `params` - params = {**params, **get_recording_channels_details(key)} - params['fs'] = params['sample_rate'] - - if acq_software == 'SpikeGLX': - spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) - spikeglx_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) - spikeglx_recording.validate_file('ap') - - if clustering_method.startswith('pykilosort'): - kilosort_triggering.run_pykilosort( - continuous_file=spikeglx_recording.root_dir / ( - spikeglx_recording.root_name + '.ap.bin'), - kilosort_output_directory=kilosort_dir, - channel_ind=params.pop('channel_ind'), - x_coords=params.pop('x_coords'), - y_coords=params.pop('y_coords'), - shank_ind=params.pop('shank_ind'), - connected=params.pop('connected'), - sample_rate=params.pop('sample_rate'), - params=params) - else: - run_kilosort = kilosort_triggering.SGLXKilosortPipeline( - npx_input_dir=spikeglx_meta_filepath.parent, - ks_output_dir=kilosort_dir, - params=params, - KS2ver=f'{Decimal(clustering_method.replace("kilosort", "")):.1f}', - run_CatGT=False) - run_kilosort.run_modules() - elif acq_software == 'Open Ephys': - oe_probe = get_openephys_probe_data(key) - - assert len(oe_probe.recording_info['recording_files']) == 1 - - # run kilosort - if clustering_method.startswith('pykilosort'): - kilosort_triggering.run_pykilosort( - continuous_file=pathlib.Path(oe_probe.recording_info['recording_files'][0]) / 'continuous.dat', - kilosort_output_directory=kilosort_dir, - channel_ind=params.pop('channel_ind'), - x_coords=params.pop('x_coords'), - y_coords=params.pop('y_coords'), - shank_ind=params.pop('shank_ind'), - connected=params.pop('connected'), - sample_rate=params.pop('sample_rate'), - params=params) - else: - run_kilosort = kilosort_triggering.OpenEphysKilosortPipeline( - npx_input_dir=oe_probe.recording_info['recording_files'][0], - ks_output_dir=kilosort_dir, - params=params, - KS2ver=f'{Decimal(clustering_method.replace("kilosort", "")):.1f}') - run_kilosort.run_modules() - else: - raise NotImplementedError(f'Automatic triggering of {clustering_method}' - f' clustering analysis is not yet supported') - - else: - raise ValueError(f'Unknown task mode: {task_mode}') - - creation_time, _, _ = kilosort.extract_clustering_info(kilosort_dir) - self.insert1({**key, 'clustering_time': creation_time}) - - -@schema -class Curation(dj.Manual): - definition = """ - # Manual curation procedure - -> Clustering - curation_id: int - --- - curation_time: datetime # time of generation of this set of curated clustering results - curation_output_dir: varchar(255) # output directory of the curated results, relative to root data directory - quality_control: bool # has this clustering result undergone quality control? - manual_curation: bool # has manual curation been performed on this clustering result? - curation_note='': varchar(2000) - """ - - def create1_from_clustering_task(self, key, curation_note=''): - """ - A function to create a new corresponding "Curation" for a particular - "ClusteringTask" - """ - if key not in Clustering(): - raise ValueError(f'No corresponding entry in Clustering available' - f' for: {key}; do `Clustering.populate(key)`') - - task_mode, output_dir = (ClusteringTask & key).fetch1( - 'task_mode', 'clustering_output_dir') - kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) - - creation_time, is_curated, is_qc = kilosort.extract_clustering_info(kilosort_dir) - # Synthesize curation_id - curation_id = dj.U().aggr(self & key, n='ifnull(max(curation_id)+1,1)').fetch1('n') - self.insert1({**key, 'curation_id': curation_id, - 'curation_time': creation_time, - 'curation_output_dir': output_dir, - 'quality_control': is_qc, - 'manual_curation': is_curated, - 'curation_note': curation_note}) - - -@schema -class CuratedClustering(dj.Imported): - definition = """ - # Clustering results of a curation. - -> Curation - """ - - class Unit(dj.Part): - definition = """ - # Properties of a given unit from a round of clustering (and curation) - -> master - unit: int - --- - -> probe.ElectrodeConfig.Electrode # electrode with highest waveform amplitude for this unit - -> ClusterQualityLabel - spike_count: int # how many spikes in this recording for this unit - spike_times: longblob # (s) spike times of this unit, relative to the start of the EphysRecording - spike_sites : longblob # array of electrode associated with each spike - spike_depths : longblob # (um) array of depths associated with each spike, relative to the (0, 0) of the probe - """ - - def make(self, key): - output_dir = (Curation & key).fetch1('curation_output_dir') - kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) - - kilosort_dataset = kilosort.Kilosort(kilosort_dir) - acq_software, sample_rate = (EphysRecording & key).fetch1( - 'acq_software', 'sampling_rate') - - sample_rate = kilosort_dataset.data['params'].get('sample_rate', sample_rate) - - # ---------- Unit ---------- - # -- Remove 0-spike units - withspike_idx = [i for i, u in enumerate(kilosort_dataset.data['cluster_ids']) - if (kilosort_dataset.data['spike_clusters'] == u).any()] - valid_units = kilosort_dataset.data['cluster_ids'][withspike_idx] - valid_unit_labels = kilosort_dataset.data['cluster_groups'][withspike_idx] - # -- Get channel and electrode-site mapping - channel2electrodes = get_neuropixels_channel2electrode_map(key, acq_software) - - # -- Spike-times -- - # spike_times_sec_adj > spike_times_sec > spike_times - spike_time_key = ('spike_times_sec_adj' if 'spike_times_sec_adj' in kilosort_dataset.data - else 'spike_times_sec' if 'spike_times_sec' - in kilosort_dataset.data else 'spike_times') - spike_times = kilosort_dataset.data[spike_time_key] - kilosort_dataset.extract_spike_depths() - - # -- Spike-sites and Spike-depths -- - spike_sites = np.array([channel2electrodes[s]['electrode'] - for s in kilosort_dataset.data['spike_sites']]) - spike_depths = kilosort_dataset.data['spike_depths'] - - # -- Insert unit, label, peak-chn - units = [] - for unit, unit_lbl in zip(valid_units, valid_unit_labels): - if (kilosort_dataset.data['spike_clusters'] == unit).any(): - unit_channel, _ = kilosort_dataset.get_best_channel(unit) - unit_spike_times = (spike_times[kilosort_dataset.data['spike_clusters'] == unit] - / sample_rate) - spike_count = len(unit_spike_times) - - units.append({ - 'unit': unit, - 'cluster_quality_label': unit_lbl, - **channel2electrodes[unit_channel], - 'spike_times': unit_spike_times, - 'spike_count': spike_count, - 'spike_sites': spike_sites[kilosort_dataset.data['spike_clusters'] == unit], - 'spike_depths': spike_depths[kilosort_dataset.data['spike_clusters'] == unit]}) - - self.insert1(key) - self.Unit.insert([{**key, **u} for u in units]) - - -@schema -class WaveformSet(dj.Imported): - definition = """ - # A set of spike waveforms for units out of a given CuratedClustering - -> CuratedClustering - """ - - class PeakWaveform(dj.Part): - definition = """ - # Mean waveform across spikes for a given unit at its representative electrode - -> master - -> CuratedClustering.Unit - --- - peak_electrode_waveform: longblob # (uV) mean waveform for a given unit at its representative electrode - """ - - class Waveform(dj.Part): - definition = """ - # Spike waveforms and their mean across spikes for the given unit - -> master - -> CuratedClustering.Unit - -> probe.ElectrodeConfig.Electrode - --- - waveform_mean: longblob # (uV) mean waveform across spikes of the given unit - waveforms=null: longblob # (uV) (spike x sample) waveforms of a sampling of spikes at the given electrode for the given unit - """ - - def make(self, key): - output_dir = (Curation & key).fetch1('curation_output_dir') - kilosort_dir = find_full_path(get_ephys_root_data_dir(), output_dir) - - kilosort_dataset = kilosort.Kilosort(kilosort_dir) - - acq_software, probe_serial_number = (EphysRecording * ProbeInsertion & key).fetch1( - 'acq_software', 'probe') - - # -- Get channel and electrode-site mapping - recording_key = (EphysRecording & key).fetch1('KEY') - channel2electrodes = get_neuropixels_channel2electrode_map(recording_key, acq_software) - - is_qc = (Curation & key).fetch1('quality_control') - - # Get all units - units = {u['unit']: u for u in (CuratedClustering.Unit & key).fetch( - as_dict=True, order_by='unit')} - - if is_qc: - unit_waveforms = np.load(kilosort_dir / 'mean_waveforms.npy') # unit x channel x sample - - def yield_unit_waveforms(): - for unit_no, unit_waveform in zip(kilosort_dataset.data['cluster_ids'], - unit_waveforms): - unit_peak_waveform = {} - unit_electrode_waveforms = [] - if unit_no in units: - for channel, channel_waveform in zip( - kilosort_dataset.data['channel_map'], - unit_waveform): - unit_electrode_waveforms.append({ - **units[unit_no], **channel2electrodes[channel], - 'waveform_mean': channel_waveform}) - if channel2electrodes[channel]['electrode'] == units[unit_no]['electrode']: - unit_peak_waveform = { - **units[unit_no], - 'peak_electrode_waveform': channel_waveform} - yield unit_peak_waveform, unit_electrode_waveforms - else: - if acq_software == 'SpikeGLX': - spikeglx_meta_filepath = get_spikeglx_meta_filepath(key) - neuropixels_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) - elif acq_software == 'Open Ephys': - session_dir = find_full_path(get_ephys_root_data_dir(), - get_session_directory(key)) - openephys_dataset = openephys.OpenEphys(session_dir) - neuropixels_recording = openephys_dataset.probes[probe_serial_number] - - def yield_unit_waveforms(): - for unit_dict in units.values(): - unit_peak_waveform = {} - unit_electrode_waveforms = [] - - spikes = unit_dict['spike_times'] - waveforms = neuropixels_recording.extract_spike_waveforms( - spikes, kilosort_dataset.data['channel_map']) # (sample x channel x spike) - waveforms = waveforms.transpose((1, 2, 0)) # (channel x spike x sample) - for channel, channel_waveform in zip( - kilosort_dataset.data['channel_map'], waveforms): - unit_electrode_waveforms.append({ - **unit_dict, **channel2electrodes[channel], - 'waveform_mean': channel_waveform.mean(axis=0), - 'waveforms': channel_waveform}) - if channel2electrodes[channel]['electrode'] == unit_dict['electrode']: - unit_peak_waveform = { - **unit_dict, - 'peak_electrode_waveform': channel_waveform.mean(axis=0)} - - yield unit_peak_waveform, unit_electrode_waveforms - - # insert waveform on a per-unit basis to mitigate potential memory issue - self.insert1(key) - for unit_peak_waveform, unit_electrode_waveforms in yield_unit_waveforms(): - if unit_peak_waveform: - self.PeakWaveform.insert1(unit_peak_waveform, ignore_extra_fields=True) - if unit_electrode_waveforms: - self.Waveform.insert(unit_electrode_waveforms, ignore_extra_fields=True) - - -# ---------------- HELPER FUNCTIONS ---------------- - -def get_spikeglx_meta_filepath(ephys_recording_key): - # attempt to retrieve from EphysRecording.EphysFile - spikeglx_meta_filepath = (EphysRecording.EphysFile & ephys_recording_key - & 'file_path LIKE "%.ap.meta"').fetch1('file_path') - - try: - spikeglx_meta_filepath = find_full_path(get_ephys_root_data_dir(), - spikeglx_meta_filepath) - except FileNotFoundError: - # if not found, search in session_dir again - if not spikeglx_meta_filepath.exists(): - session_dir = find_full_path(get_ephys_root_data_dir(), - get_session_directory( - ephys_recording_key)) - inserted_probe_serial_number = (ProbeInsertion * probe.Probe - & ephys_recording_key).fetch1('probe') - - spikeglx_meta_filepaths = [fp for fp in session_dir.rglob('*.ap.meta')] - for meta_filepath in spikeglx_meta_filepaths: - spikeglx_meta = spikeglx.SpikeGLXMeta(meta_filepath) - if str(spikeglx_meta.probe_SN) == inserted_probe_serial_number: - spikeglx_meta_filepath = meta_filepath - break - else: - raise FileNotFoundError( - 'No SpikeGLX data found for probe insertion: {}'.format(ephys_recording_key)) - - return spikeglx_meta_filepath - - -def get_openephys_probe_data(ephys_recording_key): - inserted_probe_serial_number = (ProbeInsertion * probe.Probe - & ephys_recording_key).fetch1('probe') - session_dir = find_full_path(get_ephys_root_data_dir(), - get_session_directory(ephys_recording_key)) - loaded_oe = openephys.OpenEphys(session_dir) - return loaded_oe.probes[inserted_probe_serial_number] - - -def get_neuropixels_channel2electrode_map(ephys_recording_key, acq_software): - if acq_software == 'SpikeGLX': - spikeglx_meta_filepath = get_spikeglx_meta_filepath(ephys_recording_key) - spikeglx_meta = spikeglx.SpikeGLXMeta(spikeglx_meta_filepath) - electrode_config_key = (EphysRecording * probe.ElectrodeConfig - & ephys_recording_key).fetch1('KEY') - - electrode_query = (probe.ProbeType.Electrode - * probe.ElectrodeConfig.Electrode & electrode_config_key) - - probe_electrodes = { - (shank, shank_col, shank_row): key - for key, shank, shank_col, shank_row in zip(*electrode_query.fetch( - 'KEY', 'shank', 'shank_col', 'shank_row'))} - - channel2electrode_map = { - recorded_site: probe_electrodes[(shank, shank_col, shank_row)] - for recorded_site, (shank, shank_col, shank_row, _) in enumerate( - spikeglx_meta.shankmap['data'])} - elif acq_software == 'Open Ephys': - probe_dataset = get_openephys_probe_data(ephys_recording_key) - - electrode_query = (probe.ProbeType.Electrode - * probe.ElectrodeConfig.Electrode - * EphysRecording & ephys_recording_key) - - probe_electrodes = {key['electrode']: key - for key in electrode_query.fetch('KEY')} - - channel2electrode_map = { - channel_idx: probe_electrodes[channel_idx] - for channel_idx in probe_dataset.ap_meta['channels_indices']} - - return channel2electrode_map - - -def generate_electrode_config(probe_type: str, electrodes: list): - """ - Generate and insert new ElectrodeConfig - :param probe_type: probe type (e.g. neuropixels 2.0 - SS) - :param electrodes: list of the electrode dict (keys of the probe.ProbeType.Electrode table) - :return: a dict representing a key of the probe.ElectrodeConfig table - """ - # compute hash for the electrode config (hash of dict of all ElectrodeConfig.Electrode) - electrode_config_hash = dict_to_uuid({k['electrode']: k for k in electrodes}) - - electrode_list = sorted([k['electrode'] for k in electrodes]) - electrode_gaps = ([-1] - + np.where(np.diff(electrode_list) > 1)[0].tolist() - + [len(electrode_list) - 1]) - electrode_config_name = '; '.join([ - f'{electrode_list[start + 1]}-{electrode_list[end]}' - for start, end in zip(electrode_gaps[:-1], electrode_gaps[1:])]) - - electrode_config_key = {'electrode_config_hash': electrode_config_hash} - - # ---- make new ElectrodeConfig if needed ---- - if not probe.ElectrodeConfig & electrode_config_key: - probe.ElectrodeConfig.insert1({**electrode_config_key, 'probe_type': probe_type, - 'electrode_config_name': electrode_config_name}) - probe.ElectrodeConfig.Electrode.insert({**electrode_config_key, **electrode} - for electrode in electrodes) - - return electrode_config_key - - -def get_recording_channels_details(ephys_recording_key): - channels_details = {} - - acq_software, sample_rate = (EphysRecording & ephys_recording_key).fetch1('acq_software', - 'sampling_rate') - - probe_type = (ProbeInsertion * probe.Probe & ephys_recording_key).fetch1('probe_type') - channels_details['probe_type'] = {'neuropixels 1.0 - 3A': '3A', - 'neuropixels 1.0 - 3B': 'NP1', - 'neuropixels UHD': 'NP1100', - 'neuropixels 2.0 - SS': 'NP21', - 'neuropixels 2.0 - MS': 'NP24'}[probe_type] - - electrode_config_key = (probe.ElectrodeConfig * EphysRecording & ephys_recording_key).fetch1('KEY') - channels_details['channel_ind'], channels_details['x_coords'], channels_details[ - 'y_coords'], channels_details['shank_ind'] = ( - probe.ElectrodeConfig.Electrode * probe.ProbeType.Electrode - & electrode_config_key).fetch('electrode', 'x_coord', 'y_coord', 'shank') - channels_details['sample_rate'] = sample_rate - channels_details['num_channels'] = len(channels_details['channel_ind']) - - if acq_software == 'SpikeGLX': - spikeglx_meta_filepath = get_spikeglx_meta_filepath(ephys_recording_key) - spikeglx_recording = spikeglx.SpikeGLX(spikeglx_meta_filepath.parent) - channels_details['uVPerBit'] = spikeglx_recording.get_channel_bit_volts('ap')[0] - channels_details['connected'] = np.array( - [v for *_, v in spikeglx_recording.apmeta.shankmap['data']]) - elif acq_software == 'Open Ephys': - oe_probe = get_openephys_probe_data(ephys_recording_key) - channels_details['uVPerBit'] = oe_probe.ap_meta['channels_gains'][0] - channels_details['connected'] = np.array([ - int(v == 1) for c, v in oe_probe.channels_connected.items() - if c in channels_details['channel_ind']]) - - return channels_details From 441cfe2e2ab00765ad9603d28f0bd8a50d48d1d1 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Wed, 27 Apr 2022 11:03:37 -0400 Subject: [PATCH 76/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Dimitri Yatsenko --- element_array_ephys/export/nwb/nwb.py | 7 ++----- 1 file changed, 2 insertions(+), 5 deletions(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index c0a8c8a7..6d34c6e2 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -18,11 +18,8 @@ assert probe.schema.is_activated(), 'probe not yet activated' -for ephys in (ephys_no_curation,): - if ephys.schema.is_activated(): - break -else: - raise AssertionError('ephys not yet activated') +assert ephys_no_curation.schema.is_activated, \ + "The ephys module must be activated before export." class DecimalEncoder(json.JSONEncoder): From 9ee60885e638e85746a093c791a848dbd37f2472 Mon Sep 17 00:00:00 2001 From: Ben Dichter Date: Wed, 27 Apr 2022 11:03:43 -0400 Subject: [PATCH 77/77] Update element_array_ephys/export/nwb/nwb.py Co-authored-by: Dimitri Yatsenko --- element_array_ephys/export/nwb/nwb.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/element_array_ephys/export/nwb/nwb.py b/element_array_ephys/export/nwb/nwb.py index 6d34c6e2..321159ee 100644 --- a/element_array_ephys/export/nwb/nwb.py +++ b/element_array_ephys/export/nwb/nwb.py @@ -50,7 +50,7 @@ def __init__(self, lfp_electrodes_query, chunk_length: int = 10000): self.lfp_electrodes_query & dict(electrode=self.electrodes[0]) ).fetch1(as_dict=True) - self.n_channels = len(lfp_electrodes_query) + self.n_channels = len(self.electrodes) self.n_tt = len(first_record["lfp"]) self._dtype = first_record["lfp"].dtype