From d8a8380d225564d58c350398a1a0d5db69b6348f Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Thu, 7 Dec 2023 16:27:56 -0500 Subject: [PATCH 01/53] [WIP] prototype to handle gaps --- mne/io/neuralynx/neuralynx.py | 160 +++++++++++++++++++++++++++++++--- 1 file changed, 146 insertions(+), 14 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 06d5000fcb6..774fe819e0c 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -4,13 +4,36 @@ import os import numpy as np +from numpy.testing import assert_allclose from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose +from ...annotations import Annotations from ..base import BaseRaw +from neo import AnalogSignal + +class AnalogSignalGap(object): + + def __init__(self, signal, units, sampling_rate): + + self.signal = signal + self.units = units + self.sampling_rate = sampling_rate + + def load(self, channel_indexes): + """Dummy method such that it returns object and we access .magnitude""" + + # self.magnitude = self.magnitude[channel_indexes, :] + sig = AnalogSignal(signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate) + return sig + + + @fill_doc def read_raw_neuralynx( fname, *, preload=False, exclude_fname_patterns=None, verbose=None @@ -103,25 +126,109 @@ def __init__( for segment in range(n_segments) ) - # construct an array of shape (n_total_samples,) indicating - # segment membership for each sample - sample2segment = np.concatenate( - [ - np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) - for i in range(n_segments) + # get segment start/stop times + start_times = np.array([nlx_reader.segment_t_start(block_id, i) for i in range(n_segments)]) + stop_times = np.array([nlx_reader.segment_t_stop(block_id, i) for i in range(n_segments)]) + + # find discontinuous boundaries (of length n-1) + next_start_times = start_times[1::] + previous_stop_times = stop_times[:-1] + seg_diffs = next_start_times - previous_stop_times + + # mark as discontinuous any two segments that have start/stop delta + # larger than 0.0001 sec + seg_gap_dict = {} + gap_segment_sizes = [] + delta = 1/info["sfreq"] # sampling period + gaps = ~np.isclose(seg_diffs, 0, atol=delta) + has_gaps = gaps.any() + + if has_gaps: + + logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta:.4f} \n(max = {seg_diffs[gaps].max():.3f} sec, min = {seg_diffs[gaps].min():.3f})") + + gap_starts = stop_times[:-1][gaps] # gap starts at segment offset + gap_stops = start_times[1::][gaps] # gap stops at segment onset + + # (n_gaps,) array of ints giving number of samples per inferred gap + gap_n_samps = np.array( + [len(np.arange(on_spl, off_spl)) + for on_spl, off_spl in zip(gap_starts*info["sfreq"], gap_stops*info["sfreq"]) + ] + ) + + # join + all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) + all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) + gap_indicator = np.concatenate( + [np.full(len(start_times), fill_value=0), + np.full(len(gap_starts), fill_value=1) + ] + ) + all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] + all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] + gap_indicator = gap_indicator[all_starts_ids].astype(bool) + + seg_gap_dict = { + "onsets": all_starts, # onsets in seconds + "offsets": all_stops, + "gap_n_samps": gap_n_samps, + "isgap": gap_indicator, # 0 (data segment) or 1 (invalid segment for BAD_SKIP_ACQ) + } + + gap_annotations = dict(onset=gap_starts, duration=seg_diffs[gaps], orig_time=None, description="BAD_ACQ_SKIP") + + gap_segment_sizes = [ + n for n in gap_n_samps ] + + else: + logger.info(f"All Neo segments temporally continuous at {delta:.3f} sec precision") + + # check that segment[-1] stop and segment[i] start times + # matched to microsecond precision (1e-6) + #breakpoint() + #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, + # err_msg="Segments start/end times are not temporally contiguous." + #) + + valid_segment_sizes = [ + nlx_reader.get_signal_size(block_id, i) + for i in range(n_segments) + ] + + if has_gaps: + sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[all_starts_ids] + else: + sizes_sorted = np.array(valid_segment_sizes) + + # now construct an (n_samples,) indicator variable + sample2segment = np.concatenate( + [np.full(shape=(n,), fill_value=i) + for i, n in enumerate(sizes_sorted)] ) + # construct an array of shape (n_total_samples,) indicating + # segment membership for each sample + #sample2segment = np.concatenate( + # [ + # np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) + # for i in range(n_segments) + # ] + #) + super(RawNeuralynx, self).__init__( info=info, last_samps=[n_total_samples - 1], filenames=[fname], preload=preload, - raw_extras=[dict(smp2seg=sample2segment, exclude_fnames=exclude_fnames)], + raw_extras=[dict(smp2seg=sample2segment, exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, gap_annotations=gap_annotations)], ) def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" + from quantities import Hz + from neo import Segment from neo.io import NeuralynxIO nlx_reader = NeuralynxIO( @@ -137,12 +244,14 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): ) == len(neo_block[0].segments) # collect sizes of each segment - segment_sizes = np.array( - [ - nlx_reader.get_signal_size(0, segment_id) - for segment_id in range(len(neo_block[0].segments)) - ] - ) + #segment_sizes = np.array( + # [ + # nlx_reader.get_signal_size(0, segment_id) + # for segment_id in range(len(neo_block[0].segments)) + # ] + #) + + segment_sizes = self._raw_extras[0]["segment_sizes"] # construct a (n_segments, 2) array of the first and last # sample index for each segment relative to the start of the recording @@ -188,6 +297,29 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): -1, 0 ] # express stop sample relative to segment onset + # array containing Segments + segments_arr = np.array(neo_block[0].segments, dtype=object) + + # if gaps were detected correctly insert gap Segments in between valid Segments + if self._raw_extras[0]["seg_gap_dict"]: + + gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] + gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] + + # create AnalogSignal objects representing gap data + sfreq = nlx_reader.get_signal_sampling_rate() + n_chans = np.arange(idx.start, idx.stop, idx.step).size + + for seg, n in zip(gap_segments, gap_samples): + asig = AnalogSignalGap(signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz) + seg.analogsignals.append(asig) + + n_total_segments = len(neo_block[0].segments + gap_segments) + segments_arr = np.zeros((n_total_segments,), dtype=object) + isgap = self._raw_extras[0]["seg_gap_dict"]["isgap"] + segments_arr[~isgap] = neo_block[0].segments + segments_arr[isgap] = gap_segments + # now load data from selected segments/channels via # neo.Segment.AnalogSignal.load() all_data = np.concatenate( @@ -196,7 +328,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): samples[0] : samples[-1] + 1, : ] for seg, samples in zip( - neo_block[0].segments[first_seg : last_seg + 1], sel_samples_local + segments_arr[first_seg : last_seg + 1], sel_samples_local ) for signal in seg.analogsignals ] From 67e34cf4282f09bf021080aaf9964213b3d6a91a Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:40:27 -0500 Subject: [PATCH 02/53] [ENH] define gap as > than the sampling period --- mne/io/neuralynx/neuralynx.py | 34 +++++++++++++++++++++++----------- 1 file changed, 23 insertions(+), 11 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 774fe819e0c..35eaadb8f73 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -16,7 +16,10 @@ from neo import AnalogSignal class AnalogSignalGap(object): - + """Dummy object to represent gaps in Neuralynx data as + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """ def __init__(self, signal, units, sampling_rate): self.signal = signal @@ -135,17 +138,20 @@ def __init__( previous_stop_times = stop_times[:-1] seg_diffs = next_start_times - previous_stop_times - # mark as discontinuous any two segments that have start/stop delta - # larger than 0.0001 sec + # mark as discontinuous any two segments that have + # start/stop delta larger than sampling period (1/sampling_rate) + + delta = 1/info["sfreq"] + gaps = seg_diffs > delta + has_gaps = gaps.any() + seg_gap_dict = {} gap_segment_sizes = [] - delta = 1/info["sfreq"] # sampling period - gaps = ~np.isclose(seg_diffs, 0, atol=delta) - has_gaps = gaps.any() + gap_annotations = {} if has_gaps: - logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta:.4f} \n(max = {seg_diffs[gaps].max():.3f} sec, min = {seg_diffs[gaps].min():.3f})") + logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta} sec.\n(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})") gap_starts = stop_times[:-1][gaps] # gap starts at segment offset gap_stops = start_times[1::][gaps] # gap stops at segment onset @@ -157,18 +163,23 @@ def __init__( ] ) - # join + # add the inferred gaps into the right place in the segment list all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) + + # sort the valid segment and gap times by time + all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] + all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] + + # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( [np.full(len(start_times), fill_value=0), np.full(len(gap_starts), fill_value=1) ] ) - all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] - all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] gap_indicator = gap_indicator[all_starts_ids].astype(bool) + # store this in a dict to be passed to _raw_extras seg_gap_dict = { "onsets": all_starts, # onsets in seconds "offsets": all_stops, @@ -176,6 +187,7 @@ def __init__( "isgap": gap_indicator, # 0 (data segment) or 1 (invalid segment for BAD_SKIP_ACQ) } + # TMP: annotations dict for use with mne.Annotations gap_annotations = dict(onset=gap_starts, duration=seg_diffs[gaps], orig_time=None, description="BAD_ACQ_SKIP") gap_segment_sizes = [ @@ -183,7 +195,7 @@ def __init__( ] else: - logger.info(f"All Neo segments temporally continuous at {delta:.3f} sec precision") + logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") # check that segment[-1] stop and segment[i] start times # matched to microsecond precision (1e-6) From 09c1ad94d9a28db28f6ef9173063c9600d6cefaa Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:42:13 -0500 Subject: [PATCH 03/53] [MAINT] ruff --- mne/io/neuralynx/neuralynx.py | 29 +++++++++++++---------------- 1 file changed, 13 insertions(+), 16 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 35eaadb8f73..9bb07140496 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -4,22 +4,20 @@ import os import numpy as np -from numpy.testing import assert_allclose +from neo import AnalogSignal from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose -from ...annotations import Annotations from ..base import BaseRaw -from neo import AnalogSignal - class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. """ + def __init__(self, signal, units, sampling_rate): self.signal = signal @@ -28,11 +26,10 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" - # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], + sig = AnalogSignal(signal=self.signal[channel_indexes, :], units=self.units, - sampling_rate=self.sampling_rate) + sampling_rate=self.sampling_rate) return sig @@ -138,10 +135,10 @@ def __init__( previous_stop_times = stop_times[:-1] seg_diffs = next_start_times - previous_stop_times - # mark as discontinuous any two segments that have + # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) - delta = 1/info["sfreq"] + delta = 1/info["sfreq"] gaps = seg_diffs > delta has_gaps = gaps.any() @@ -158,7 +155,7 @@ def __init__( # (n_gaps,) array of ints giving number of samples per inferred gap gap_n_samps = np.array( - [len(np.arange(on_spl, off_spl)) + [len(np.arange(on_spl, off_spl)) for on_spl, off_spl in zip(gap_starts*info["sfreq"], gap_stops*info["sfreq"]) ] ) @@ -166,11 +163,11 @@ def __init__( # add the inferred gaps into the right place in the segment list all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) - + # sort the valid segment and gap times by time all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] - + # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( [np.full(len(start_times), fill_value=0), @@ -208,7 +205,7 @@ def __init__( nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] - + if has_gaps: sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[all_starts_ids] else: @@ -239,9 +236,9 @@ def __init__( def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" - from quantities import Hz from neo import Segment from neo.io import NeuralynxIO + from quantities import Hz nlx_reader = NeuralynxIO( dirname=self._filenames[fi], From e9c2de2d6bca7011443ccc004d0c0508a21e0bfa Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:19:53 -0500 Subject: [PATCH 04/53] [FIX] ruff formatting pt2 --- mne/io/neuralynx/neuralynx.py | 97 ++++++++++++++++++++++------------- 1 file changed, 61 insertions(+), 36 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 9bb07140496..4d267f21e75 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -19,7 +19,6 @@ class AnalogSignalGap(object): """ def __init__(self, signal, units, sampling_rate): - self.signal = signal self.units = units self.sampling_rate = sampling_rate @@ -27,13 +26,14 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], - units=self.units, - sampling_rate=self.sampling_rate) + sig = AnalogSignal( + signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate, + ) return sig - @fill_doc def read_raw_neuralynx( fname, *, preload=False, exclude_fname_patterns=None, verbose=None @@ -127,8 +127,12 @@ def __init__( ) # get segment start/stop times - start_times = np.array([nlx_reader.segment_t_start(block_id, i) for i in range(n_segments)]) - stop_times = np.array([nlx_reader.segment_t_stop(block_id, i) for i in range(n_segments)]) + start_times = np.array( + [nlx_reader.segment_t_start(block_id, i) for i in range(n_segments)] + ) + stop_times = np.array( + [nlx_reader.segment_t_stop(block_id, i) for i in range(n_segments)] + ) # find discontinuous boundaries (of length n-1) next_start_times = start_times[1::] @@ -138,7 +142,7 @@ def __init__( # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) - delta = 1/info["sfreq"] + delta = 1 / info["sfreq"] gaps = seg_diffs > delta has_gaps = gaps.any() @@ -147,16 +151,22 @@ def __init__( gap_annotations = {} if has_gaps: - - logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta} sec.\n(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})") + logger.info( + f"N = {gaps.sum()} discontinuous Neo segments detected " + + "with delta > {delta} sec. +\n" + + "(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" + ) gap_starts = stop_times[:-1][gaps] # gap starts at segment offset gap_stops = start_times[1::][gaps] # gap stops at segment onset # (n_gaps,) array of ints giving number of samples per inferred gap gap_n_samps = np.array( - [len(np.arange(on_spl, off_spl)) - for on_spl, off_spl in zip(gap_starts*info["sfreq"], gap_stops*info["sfreq"]) + [ + len(np.arange(on_spl, off_spl)) + for on_spl, off_spl in zip( + gap_starts * info["sfreq"], gap_stops * info["sfreq"] + ) ] ) @@ -170,68 +180,82 @@ def __init__( # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( - [np.full(len(start_times), fill_value=0), - np.full(len(gap_starts), fill_value=1) + [ + np.full(len(start_times), fill_value=0), + np.full(len(gap_starts), fill_value=1), ] ) gap_indicator = gap_indicator[all_starts_ids].astype(bool) # store this in a dict to be passed to _raw_extras seg_gap_dict = { - "onsets": all_starts, # onsets in seconds + "onsets": all_starts, # onsets in seconds "offsets": all_stops, "gap_n_samps": gap_n_samps, - "isgap": gap_indicator, # 0 (data segment) or 1 (invalid segment for BAD_SKIP_ACQ) + "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) } # TMP: annotations dict for use with mne.Annotations - gap_annotations = dict(onset=gap_starts, duration=seg_diffs[gaps], orig_time=None, description="BAD_ACQ_SKIP") + gap_annotations = dict( + onset=gap_starts, + duration=seg_diffs[gaps], + orig_time=None, + description="BAD_ACQ_SKIP", + ) - gap_segment_sizes = [ - n for n in gap_n_samps - ] + gap_segment_sizes = [n for n in gap_n_samps] else: - logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") + logger.info( + f"All Neo segments temporally continuous at {delta} sec precision." + ) # check that segment[-1] stop and segment[i] start times # matched to microsecond precision (1e-6) - #breakpoint() - #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, + # breakpoint() + # assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, # err_msg="Segments start/end times are not temporally contiguous." - #) + # ) valid_segment_sizes = [ - nlx_reader.get_signal_size(block_id, i) - for i in range(n_segments) + nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] if has_gaps: - sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[all_starts_ids] + sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[ + all_starts_ids + ] else: sizes_sorted = np.array(valid_segment_sizes) # now construct an (n_samples,) indicator variable sample2segment = np.concatenate( - [np.full(shape=(n,), fill_value=i) - for i, n in enumerate(sizes_sorted)] + [np.full(shape=(n,), fill_value=i) for i, n in enumerate(sizes_sorted)] ) # construct an array of shape (n_total_samples,) indicating # segment membership for each sample - #sample2segment = np.concatenate( + # sample2segment = np.concatenate( # [ # np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) # for i in range(n_segments) # ] - #) + # ) super(RawNeuralynx, self).__init__( info=info, last_samps=[n_total_samples - 1], filenames=[fname], preload=preload, - raw_extras=[dict(smp2seg=sample2segment, exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, gap_annotations=gap_annotations)], + raw_extras=[ + dict( + smp2seg=sample2segment, + exclude_fnames=exclude_fnames, + segment_sizes=sizes_sorted, + seg_gap_dict=seg_gap_dict, + gap_annotations=gap_annotations, + ) + ], ) def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): @@ -253,12 +277,12 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): ) == len(neo_block[0].segments) # collect sizes of each segment - #segment_sizes = np.array( + # segment_sizes = np.array( # [ # nlx_reader.get_signal_size(0, segment_id) # for segment_id in range(len(neo_block[0].segments)) # ] - #) + # ) segment_sizes = self._raw_extras[0]["segment_sizes"] @@ -311,7 +335,6 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): # if gaps were detected correctly insert gap Segments in between valid Segments if self._raw_extras[0]["seg_gap_dict"]: - gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] @@ -320,7 +343,9 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): n_chans = np.arange(idx.start, idx.stop, idx.step).size for seg, n in zip(gap_segments, gap_samples): - asig = AnalogSignalGap(signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz) + asig = AnalogSignalGap( + signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz + ) seg.analogsignals.append(asig) n_total_segments = len(neo_block[0].segments + gap_segments) From bf1cfa0a0d65aa8d1516873346a64fca4c8b6fb8 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:26:41 -0500 Subject: [PATCH 05/53] [FIX] more style fixes --- mne/io/neuralynx/neuralynx.py | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 4d267f21e75..5022f18be46 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -13,9 +13,11 @@ class AnalogSignalGap(object): - """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """Dummy object to represent gaps in Neuralynx data. + + Creates a AnalogSignalProxy-like object. + Propagate `signal`, `units`, and `sampling_rate` attributes + to the `AnalogSignal` object returned by `load()`. """ def __init__(self, signal, units, sampling_rate): @@ -24,7 +26,7 @@ def __init__(self, signal, units, sampling_rate): self.sampling_rate = sampling_rate def load(self, channel_indexes): - """Dummy method such that it returns object and we access .magnitude""" + """Return AnalogSignal object.""" # self.magnitude = self.magnitude[channel_indexes, :] sig = AnalogSignal( signal=self.signal[channel_indexes, :], @@ -210,13 +212,6 @@ def __init__( f"All Neo segments temporally continuous at {delta} sec precision." ) - # check that segment[-1] stop and segment[i] start times - # matched to microsecond precision (1e-6) - # breakpoint() - # assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, - # err_msg="Segments start/end times are not temporally contiguous." - # ) - valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] From 7fc96f1381178573ae252121cc97e82b636473a7 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Thu, 7 Dec 2023 16:27:56 -0500 Subject: [PATCH 06/53] [WIP] prototype to handle gaps --- mne/io/neuralynx/neuralynx.py | 22 ++++++++++++++++++++++ 1 file changed, 22 insertions(+) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 5022f18be46..4a062b8e42b 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -9,9 +9,31 @@ from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose +from ...annotations import Annotations from ..base import BaseRaw +from neo import AnalogSignal + +class AnalogSignalGap(object): + + def __init__(self, signal, units, sampling_rate): + + self.signal = signal + self.units = units + self.sampling_rate = sampling_rate + + def load(self, channel_indexes): + """Dummy method such that it returns object and we access .magnitude""" + + # self.magnitude = self.magnitude[channel_indexes, :] + sig = AnalogSignal(signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate) + return sig + + + class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data. From d57c56df869a7ae7b135def59fa46b034a71f8bf Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:40:27 -0500 Subject: [PATCH 07/53] [ENH] define gap as > than the sampling period --- mne/io/neuralynx/neuralynx.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 4a062b8e42b..686c555f4f4 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -16,7 +16,10 @@ from neo import AnalogSignal class AnalogSignalGap(object): - + """Dummy object to represent gaps in Neuralynx data as + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """ def __init__(self, signal, units, sampling_rate): self.signal = signal From cd4a69a459ff1c2a18a55aae506c93828f777527 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:42:13 -0500 Subject: [PATCH 08/53] [MAINT] ruff --- mne/io/neuralynx/neuralynx.py | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 686c555f4f4..54b4662de34 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -9,17 +9,15 @@ from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose -from ...annotations import Annotations from ..base import BaseRaw -from neo import AnalogSignal - class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. """ + def __init__(self, signal, units, sampling_rate): self.signal = signal @@ -28,11 +26,10 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" - # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], + sig = AnalogSignal(signal=self.signal[channel_indexes, :], units=self.units, - sampling_rate=self.sampling_rate) + sampling_rate=self.sampling_rate) return sig From 07d64eb322b70f6e5868520113edbdcf3a51ea39 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:19:53 -0500 Subject: [PATCH 09/53] [FIX] ruff formatting pt2 --- mne/io/neuralynx/neuralynx.py | 20 +++++++++++++------- 1 file changed, 13 insertions(+), 7 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 54b4662de34..277c2759836 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -19,7 +19,6 @@ class AnalogSignalGap(object): """ def __init__(self, signal, units, sampling_rate): - self.signal = signal self.units = units self.sampling_rate = sampling_rate @@ -27,9 +26,11 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], - units=self.units, - sampling_rate=self.sampling_rate) + sig = AnalogSignal( + signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate, + ) return sig @@ -230,9 +231,14 @@ def __init__( gap_segment_sizes = [n for n in gap_n_samps] else: - logger.info( - f"All Neo segments temporally continuous at {delta} sec precision." - ) + logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") + + # check that segment[-1] stop and segment[i] start times + # matched to microsecond precision (1e-6) + #breakpoint() + #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, + # err_msg="Segments start/end times are not temporally contiguous." + #) valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) From 4357ab366dce824ee9305420d39d835217d33b9b Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:26:41 -0500 Subject: [PATCH 10/53] [FIX] more style fixes --- mne/io/neuralynx/neuralynx.py | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 277c2759836..85943ec0e14 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -13,9 +13,11 @@ class AnalogSignalGap(object): - """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """Dummy object to represent gaps in Neuralynx data. + + Creates a AnalogSignalProxy-like object. + Propagate `signal`, `units`, and `sampling_rate` attributes + to the `AnalogSignal` object returned by `load()`. """ def __init__(self, signal, units, sampling_rate): @@ -24,7 +26,7 @@ def __init__(self, signal, units, sampling_rate): self.sampling_rate = sampling_rate def load(self, channel_indexes): - """Dummy method such that it returns object and we access .magnitude""" + """Return AnalogSignal object.""" # self.magnitude = self.magnitude[channel_indexes, :] sig = AnalogSignal( signal=self.signal[channel_indexes, :], @@ -233,13 +235,6 @@ def __init__( else: logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") - # check that segment[-1] stop and segment[i] start times - # matched to microsecond precision (1e-6) - #breakpoint() - #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, - # err_msg="Segments start/end times are not temporally contiguous." - #) - valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] From aedc9a4849688e4e1680b9b10fd0e10eeb0f0112 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 11 Dec 2023 16:21:53 -0500 Subject: [PATCH 11/53] [FIX] add annotations dict to _raw_extras --- mne/io/neuralynx/neuralynx.py | 82 +++++++++++++++++------------------ 1 file changed, 40 insertions(+), 42 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 85943ec0e14..a6af9fe9442 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -18,6 +18,12 @@ class AnalogSignalGap(object): Creates a AnalogSignalProxy-like object. Propagate `signal`, `units`, and `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + + Parameters + ---------- + signal : array-like + Array of shape (n_channels, n_samples) containing the data. + """ def __init__(self, signal, units, sampling_rate): @@ -29,7 +35,7 @@ def load(self, channel_indexes): """Return AnalogSignal object.""" # self.magnitude = self.magnitude[channel_indexes, :] sig = AnalogSignal( - signal=self.signal[channel_indexes, :], + signal=self.signal[:, channel_indexes], units=self.units, sampling_rate=self.sampling_rate, ) @@ -148,10 +154,6 @@ def __init__( # the sample sizes of all segments n_segments = nlx_reader.header["nb_segment"][0] block_id = 0 # assumes there's only one block of recording - n_total_samples = sum( - nlx_reader.get_signal_size(block_id, segment) - for segment in range(n_segments) - ) # get segment start/stop times start_times = np.array( @@ -168,20 +170,19 @@ def __init__( # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) - + logger.info("Checking for temporal discontinuities in Neo data segments.") delta = 1 / info["sfreq"] gaps = seg_diffs > delta has_gaps = gaps.any() seg_gap_dict = {} - gap_segment_sizes = [] - gap_annotations = {} if has_gaps: logger.info( - f"N = {gaps.sum()} discontinuous Neo segments detected " + - "with delta > {delta} sec. +\n" + - "(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" + f"N = {gaps.sum()} discontinuous Neo segments detected " + + f"with delta > {delta} sec.\n" + + f"(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" + + "\nAnnotating gaps as BAD_ACQ_SKIP." ) gap_starts = stop_times[:-1][gaps] # gap starts at segment offset @@ -197,13 +198,8 @@ def __init__( ] ) - # add the inferred gaps into the right place in the segment list + # get sort indices for all segments (valid and gap) in ascending order all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) - all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) - - # sort the valid segment and gap times by time - all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] - all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( @@ -216,22 +212,10 @@ def __init__( # store this in a dict to be passed to _raw_extras seg_gap_dict = { - "onsets": all_starts, # onsets in seconds - "offsets": all_stops, "gap_n_samps": gap_n_samps, "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) } - # TMP: annotations dict for use with mne.Annotations - gap_annotations = dict( - onset=gap_starts, - duration=seg_diffs[gaps], - orig_time=None, - description="BAD_ACQ_SKIP", - ) - - gap_segment_sizes = [n for n in gap_n_samps] - else: logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") @@ -240,7 +224,7 @@ def __init__( ] if has_gaps: - sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[ + sizes_sorted = np.concatenate([valid_segment_sizes, gap_n_samps])[ all_starts_ids ] else: @@ -251,18 +235,30 @@ def __init__( [np.full(shape=(n,), fill_value=i) for i, n in enumerate(sizes_sorted)] ) - # construct an array of shape (n_total_samples,) indicating - # segment membership for each sample - # sample2segment = np.concatenate( - # [ - # np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) - # for i in range(n_segments) - # ] - # ) + if has_gaps: + gap_seg_ids = np.unique(sample2segment)[gap_indicator] + gap_start_ids = np.array( + [np.where(sample2segment == seg_id)[0][0] for seg_id in gap_seg_ids] + ) + + # recreate time axis for gap annotations + mne_times = np.arange(0, len(sample2segment)) / info["sfreq"] + + assert len(gap_start_ids) == len(gap_n_samps) + annotations = dict( + onset=[mne_times[onset_id] for onset_id in gap_start_ids], + duration=[ + mne_times[onset_id + (n - 1)] - mne_times[onset_id] + for onset_id, n in zip(gap_start_ids, gap_n_samps) + ], + descriptions=["BAD_ACQ_SKIP" for _ in gap_start_ids], + ) + else: + annotations = None super(RawNeuralynx, self).__init__( info=info, - last_samps=[n_total_samples - 1], + last_samps=[sizes_sorted.sum() - 1], filenames=[fname], preload=preload, raw_extras=[ @@ -271,7 +267,7 @@ def __init__( exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, - gap_annotations=gap_annotations, + gap_annotations=annotations, ) ], ) @@ -356,7 +352,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] - # create AnalogSignal objects representing gap data + # create AnalogSignal objects representing gap data filled with 0's sfreq = nlx_reader.get_signal_sampling_rate() n_chans = np.arange(idx.start, idx.stop, idx.step).size @@ -368,12 +364,14 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): n_total_segments = len(neo_block[0].segments + gap_segments) segments_arr = np.zeros((n_total_segments,), dtype=object) + + # insert inferred gap segments at the right place in between valid segments isgap = self._raw_extras[0]["seg_gap_dict"]["isgap"] segments_arr[~isgap] = neo_block[0].segments segments_arr[isgap] = gap_segments # now load data from selected segments/channels via - # neo.Segment.AnalogSignal.load() + # neo.Segment.AnalogSignal.load() or AnalogSignalGap.load() all_data = np.concatenate( [ signal.load(channel_indexes=idx).magnitude[ From 76a709c290db1d2e99471c6c3336c58e8a394d0b Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 11 Dec 2023 16:46:26 -0500 Subject: [PATCH 12/53] [MAINT] remove commented section --- mne/io/neuralynx/neuralynx.py | 8 -------- 1 file changed, 8 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index a6af9fe9442..ab24ecc2dd7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -290,14 +290,6 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): [len(segment.analogsignals) for segment in neo_block[0].segments] ) == len(neo_block[0].segments) - # collect sizes of each segment - # segment_sizes = np.array( - # [ - # nlx_reader.get_signal_size(0, segment_id) - # for segment_id in range(len(neo_block[0].segments)) - # ] - # ) - segment_sizes = self._raw_extras[0]["segment_sizes"] # construct a (n_segments, 2) array of the first and last From 593981e17a7c5734cdc7532e62a0dbdae976f904 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 08:58:24 -0500 Subject: [PATCH 13/53] [FIX] remove doubled code from conflicts --- mne/io/neuralynx/neuralynx.py | 29 +++-------------------------- 1 file changed, 3 insertions(+), 26 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index ab24ecc2dd7..476a319a6b7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -42,31 +42,6 @@ def load(self, channel_indexes): return sig - -class AnalogSignalGap(object): - """Dummy object to represent gaps in Neuralynx data. - - Creates a AnalogSignalProxy-like object. - Propagate `signal`, `units`, and `sampling_rate` attributes - to the `AnalogSignal` object returned by `load()`. - """ - - def __init__(self, signal, units, sampling_rate): - self.signal = signal - self.units = units - self.sampling_rate = sampling_rate - - def load(self, channel_indexes): - """Return AnalogSignal object.""" - # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal( - signal=self.signal[channel_indexes, :], - units=self.units, - sampling_rate=self.sampling_rate, - ) - return sig - - @fill_doc def read_raw_neuralynx( fname, *, preload=False, exclude_fname_patterns=None, verbose=None @@ -217,7 +192,9 @@ def __init__( } else: - logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") + logger.info( + f"All Neo segments temporally continuous at {delta} sec precision." + ) valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) From 550dab8619981eb7d3e56a52f704d16740649546 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 09:02:12 -0500 Subject: [PATCH 14/53] [MAINT] fill docstring --- mne/io/neuralynx/neuralynx.py | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 476a319a6b7..81cc3f6f80a 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -17,13 +17,21 @@ class AnalogSignalGap(object): Creates a AnalogSignalProxy-like object. Propagate `signal`, `units`, and `sampling_rate` attributes - to the `AnalogSignal` object returned by `load()`. + to the `AnalogSignal` init returned by `load()`. Parameters ---------- signal : array-like Array of shape (n_channels, n_samples) containing the data. + units : str + Units of the data. (e.g., 'uV') + sampling_rate : quantity + Sampling rate of the data. (e.g., 4000 * pq.Hz) + Returns + ------- + sig : instance of AnalogSignal + A AnalogSignal object representing a gap in Neuralynx data. """ def __init__(self, signal, units, sampling_rate): From 98ea12926bf6711e366d8e0fa1076d4f18930400 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 09:45:58 -0500 Subject: [PATCH 15/53] [FIX] move neo imports at the top --- mne/io/neuralynx/neuralynx.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 81cc3f6f80a..c825f8195f5 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -4,13 +4,16 @@ import os import numpy as np -from neo import AnalogSignal from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose from ..base import BaseRaw +_soft_import("neo", "Reading NeuralynxIO files", strict=True) +from neo import AnalogSignal +from neo.io import NeuralynxIO + class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data. @@ -98,9 +101,6 @@ def __init__( exclude_fname_patterns=None, verbose=None, ): - _soft_import("neo", "Reading NeuralynxIO files", strict=True) - from neo.io import NeuralynxIO - fname = _check_fname(fname, "read", True, "fname", need_dir=True) logger.info(f"Checking files in {fname}") From 135576306faf7eb9062937d254836225a380449d Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 09:51:36 -0500 Subject: [PATCH 16/53] [FIX] move neo soft imports within methods --- mne/io/neuralynx/neuralynx.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index c825f8195f5..2d51a13ed08 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -10,10 +10,6 @@ from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose from ..base import BaseRaw -_soft_import("neo", "Reading NeuralynxIO files", strict=True) -from neo import AnalogSignal -from neo.io import NeuralynxIO - class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data. @@ -44,7 +40,9 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Return AnalogSignal object.""" - # self.magnitude = self.magnitude[channel_indexes, :] + _soft_import("neo", "Reading NeuralynxIO files", strict=True) + from neo import AnalogSignal + sig = AnalogSignal( signal=self.signal[:, channel_indexes], units=self.units, @@ -103,6 +101,9 @@ def __init__( ): fname = _check_fname(fname, "read", True, "fname", need_dir=True) + _soft_import("neo", "Reading NeuralynxIO files", strict=True) + from neo.io import NeuralynxIO + logger.info(f"Checking files in {fname}") # construct a list of filenames to ignore From a32ff97c8e0b7b38910cf8928a44854c3695ea06 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 16:30:22 -0500 Subject: [PATCH 17/53] [ENH] set annotations post init --- mne/io/neuralynx/neuralynx.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 2d51a13ed08..7da2ccda089 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -7,6 +7,7 @@ from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one +from ...annotations import Annotations from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose from ..base import BaseRaw @@ -221,6 +222,7 @@ def __init__( [np.full(shape=(n,), fill_value=i) for i, n in enumerate(sizes_sorted)] ) + # construct Annotations() if has_gaps: gap_seg_ids = np.unique(sample2segment)[gap_indicator] gap_start_ids = np.array( @@ -231,13 +233,13 @@ def __init__( mne_times = np.arange(0, len(sample2segment)) / info["sfreq"] assert len(gap_start_ids) == len(gap_n_samps) - annotations = dict( + annotations = Annotations( onset=[mne_times[onset_id] for onset_id in gap_start_ids], duration=[ mne_times[onset_id + (n - 1)] - mne_times[onset_id] for onset_id, n in zip(gap_start_ids, gap_n_samps) ], - descriptions=["BAD_ACQ_SKIP" for _ in gap_start_ids], + description=["BAD_ACQ_SKIP" for _ in gap_start_ids], ) else: annotations = None @@ -253,11 +255,12 @@ def __init__( exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, - gap_annotations=annotations, ) ], ) + self.set_annotations(annotations) + def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" from neo import Segment From b633be41e007c20b11a077a6758df52b7e49e2c1 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Thu, 7 Dec 2023 16:27:56 -0500 Subject: [PATCH 18/53] [WIP] prototype to handle gaps --- mne/io/neuralynx/neuralynx.py | 160 +++++++++++++++++++++++++++++++--- 1 file changed, 146 insertions(+), 14 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 06d5000fcb6..774fe819e0c 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -4,13 +4,36 @@ import os import numpy as np +from numpy.testing import assert_allclose from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose +from ...annotations import Annotations from ..base import BaseRaw +from neo import AnalogSignal + +class AnalogSignalGap(object): + + def __init__(self, signal, units, sampling_rate): + + self.signal = signal + self.units = units + self.sampling_rate = sampling_rate + + def load(self, channel_indexes): + """Dummy method such that it returns object and we access .magnitude""" + + # self.magnitude = self.magnitude[channel_indexes, :] + sig = AnalogSignal(signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate) + return sig + + + @fill_doc def read_raw_neuralynx( fname, *, preload=False, exclude_fname_patterns=None, verbose=None @@ -103,25 +126,109 @@ def __init__( for segment in range(n_segments) ) - # construct an array of shape (n_total_samples,) indicating - # segment membership for each sample - sample2segment = np.concatenate( - [ - np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) - for i in range(n_segments) + # get segment start/stop times + start_times = np.array([nlx_reader.segment_t_start(block_id, i) for i in range(n_segments)]) + stop_times = np.array([nlx_reader.segment_t_stop(block_id, i) for i in range(n_segments)]) + + # find discontinuous boundaries (of length n-1) + next_start_times = start_times[1::] + previous_stop_times = stop_times[:-1] + seg_diffs = next_start_times - previous_stop_times + + # mark as discontinuous any two segments that have start/stop delta + # larger than 0.0001 sec + seg_gap_dict = {} + gap_segment_sizes = [] + delta = 1/info["sfreq"] # sampling period + gaps = ~np.isclose(seg_diffs, 0, atol=delta) + has_gaps = gaps.any() + + if has_gaps: + + logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta:.4f} \n(max = {seg_diffs[gaps].max():.3f} sec, min = {seg_diffs[gaps].min():.3f})") + + gap_starts = stop_times[:-1][gaps] # gap starts at segment offset + gap_stops = start_times[1::][gaps] # gap stops at segment onset + + # (n_gaps,) array of ints giving number of samples per inferred gap + gap_n_samps = np.array( + [len(np.arange(on_spl, off_spl)) + for on_spl, off_spl in zip(gap_starts*info["sfreq"], gap_stops*info["sfreq"]) + ] + ) + + # join + all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) + all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) + gap_indicator = np.concatenate( + [np.full(len(start_times), fill_value=0), + np.full(len(gap_starts), fill_value=1) + ] + ) + all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] + all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] + gap_indicator = gap_indicator[all_starts_ids].astype(bool) + + seg_gap_dict = { + "onsets": all_starts, # onsets in seconds + "offsets": all_stops, + "gap_n_samps": gap_n_samps, + "isgap": gap_indicator, # 0 (data segment) or 1 (invalid segment for BAD_SKIP_ACQ) + } + + gap_annotations = dict(onset=gap_starts, duration=seg_diffs[gaps], orig_time=None, description="BAD_ACQ_SKIP") + + gap_segment_sizes = [ + n for n in gap_n_samps ] + + else: + logger.info(f"All Neo segments temporally continuous at {delta:.3f} sec precision") + + # check that segment[-1] stop and segment[i] start times + # matched to microsecond precision (1e-6) + #breakpoint() + #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, + # err_msg="Segments start/end times are not temporally contiguous." + #) + + valid_segment_sizes = [ + nlx_reader.get_signal_size(block_id, i) + for i in range(n_segments) + ] + + if has_gaps: + sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[all_starts_ids] + else: + sizes_sorted = np.array(valid_segment_sizes) + + # now construct an (n_samples,) indicator variable + sample2segment = np.concatenate( + [np.full(shape=(n,), fill_value=i) + for i, n in enumerate(sizes_sorted)] ) + # construct an array of shape (n_total_samples,) indicating + # segment membership for each sample + #sample2segment = np.concatenate( + # [ + # np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) + # for i in range(n_segments) + # ] + #) + super(RawNeuralynx, self).__init__( info=info, last_samps=[n_total_samples - 1], filenames=[fname], preload=preload, - raw_extras=[dict(smp2seg=sample2segment, exclude_fnames=exclude_fnames)], + raw_extras=[dict(smp2seg=sample2segment, exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, gap_annotations=gap_annotations)], ) def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" + from quantities import Hz + from neo import Segment from neo.io import NeuralynxIO nlx_reader = NeuralynxIO( @@ -137,12 +244,14 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): ) == len(neo_block[0].segments) # collect sizes of each segment - segment_sizes = np.array( - [ - nlx_reader.get_signal_size(0, segment_id) - for segment_id in range(len(neo_block[0].segments)) - ] - ) + #segment_sizes = np.array( + # [ + # nlx_reader.get_signal_size(0, segment_id) + # for segment_id in range(len(neo_block[0].segments)) + # ] + #) + + segment_sizes = self._raw_extras[0]["segment_sizes"] # construct a (n_segments, 2) array of the first and last # sample index for each segment relative to the start of the recording @@ -188,6 +297,29 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): -1, 0 ] # express stop sample relative to segment onset + # array containing Segments + segments_arr = np.array(neo_block[0].segments, dtype=object) + + # if gaps were detected correctly insert gap Segments in between valid Segments + if self._raw_extras[0]["seg_gap_dict"]: + + gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] + gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] + + # create AnalogSignal objects representing gap data + sfreq = nlx_reader.get_signal_sampling_rate() + n_chans = np.arange(idx.start, idx.stop, idx.step).size + + for seg, n in zip(gap_segments, gap_samples): + asig = AnalogSignalGap(signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz) + seg.analogsignals.append(asig) + + n_total_segments = len(neo_block[0].segments + gap_segments) + segments_arr = np.zeros((n_total_segments,), dtype=object) + isgap = self._raw_extras[0]["seg_gap_dict"]["isgap"] + segments_arr[~isgap] = neo_block[0].segments + segments_arr[isgap] = gap_segments + # now load data from selected segments/channels via # neo.Segment.AnalogSignal.load() all_data = np.concatenate( @@ -196,7 +328,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): samples[0] : samples[-1] + 1, : ] for seg, samples in zip( - neo_block[0].segments[first_seg : last_seg + 1], sel_samples_local + segments_arr[first_seg : last_seg + 1], sel_samples_local ) for signal in seg.analogsignals ] From 350edbc59e0228451fbdaab1afe6fc066db0fff9 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:40:27 -0500 Subject: [PATCH 19/53] [ENH] define gap as > than the sampling period --- mne/io/neuralynx/neuralynx.py | 34 +++++++++++++++++++++++----------- 1 file changed, 23 insertions(+), 11 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 774fe819e0c..35eaadb8f73 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -16,7 +16,10 @@ from neo import AnalogSignal class AnalogSignalGap(object): - + """Dummy object to represent gaps in Neuralynx data as + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """ def __init__(self, signal, units, sampling_rate): self.signal = signal @@ -135,17 +138,20 @@ def __init__( previous_stop_times = stop_times[:-1] seg_diffs = next_start_times - previous_stop_times - # mark as discontinuous any two segments that have start/stop delta - # larger than 0.0001 sec + # mark as discontinuous any two segments that have + # start/stop delta larger than sampling period (1/sampling_rate) + + delta = 1/info["sfreq"] + gaps = seg_diffs > delta + has_gaps = gaps.any() + seg_gap_dict = {} gap_segment_sizes = [] - delta = 1/info["sfreq"] # sampling period - gaps = ~np.isclose(seg_diffs, 0, atol=delta) - has_gaps = gaps.any() + gap_annotations = {} if has_gaps: - logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta:.4f} \n(max = {seg_diffs[gaps].max():.3f} sec, min = {seg_diffs[gaps].min():.3f})") + logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta} sec.\n(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})") gap_starts = stop_times[:-1][gaps] # gap starts at segment offset gap_stops = start_times[1::][gaps] # gap stops at segment onset @@ -157,18 +163,23 @@ def __init__( ] ) - # join + # add the inferred gaps into the right place in the segment list all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) + + # sort the valid segment and gap times by time + all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] + all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] + + # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( [np.full(len(start_times), fill_value=0), np.full(len(gap_starts), fill_value=1) ] ) - all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] - all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] gap_indicator = gap_indicator[all_starts_ids].astype(bool) + # store this in a dict to be passed to _raw_extras seg_gap_dict = { "onsets": all_starts, # onsets in seconds "offsets": all_stops, @@ -176,6 +187,7 @@ def __init__( "isgap": gap_indicator, # 0 (data segment) or 1 (invalid segment for BAD_SKIP_ACQ) } + # TMP: annotations dict for use with mne.Annotations gap_annotations = dict(onset=gap_starts, duration=seg_diffs[gaps], orig_time=None, description="BAD_ACQ_SKIP") gap_segment_sizes = [ @@ -183,7 +195,7 @@ def __init__( ] else: - logger.info(f"All Neo segments temporally continuous at {delta:.3f} sec precision") + logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") # check that segment[-1] stop and segment[i] start times # matched to microsecond precision (1e-6) From 38e457b4e276282e9cf670e771635983e5df8e27 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:42:13 -0500 Subject: [PATCH 20/53] [MAINT] ruff --- mne/io/neuralynx/neuralynx.py | 29 +++++++++++++---------------- 1 file changed, 13 insertions(+), 16 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 35eaadb8f73..9bb07140496 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -4,22 +4,20 @@ import os import numpy as np -from numpy.testing import assert_allclose +from neo import AnalogSignal from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose -from ...annotations import Annotations from ..base import BaseRaw -from neo import AnalogSignal - class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. """ + def __init__(self, signal, units, sampling_rate): self.signal = signal @@ -28,11 +26,10 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" - # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], + sig = AnalogSignal(signal=self.signal[channel_indexes, :], units=self.units, - sampling_rate=self.sampling_rate) + sampling_rate=self.sampling_rate) return sig @@ -138,10 +135,10 @@ def __init__( previous_stop_times = stop_times[:-1] seg_diffs = next_start_times - previous_stop_times - # mark as discontinuous any two segments that have + # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) - delta = 1/info["sfreq"] + delta = 1/info["sfreq"] gaps = seg_diffs > delta has_gaps = gaps.any() @@ -158,7 +155,7 @@ def __init__( # (n_gaps,) array of ints giving number of samples per inferred gap gap_n_samps = np.array( - [len(np.arange(on_spl, off_spl)) + [len(np.arange(on_spl, off_spl)) for on_spl, off_spl in zip(gap_starts*info["sfreq"], gap_stops*info["sfreq"]) ] ) @@ -166,11 +163,11 @@ def __init__( # add the inferred gaps into the right place in the segment list all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) - + # sort the valid segment and gap times by time all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] - + # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( [np.full(len(start_times), fill_value=0), @@ -208,7 +205,7 @@ def __init__( nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] - + if has_gaps: sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[all_starts_ids] else: @@ -239,9 +236,9 @@ def __init__( def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" - from quantities import Hz from neo import Segment from neo.io import NeuralynxIO + from quantities import Hz nlx_reader = NeuralynxIO( dirname=self._filenames[fi], From 647647f44d9e2827f14fec504d1973846e193cae Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:19:53 -0500 Subject: [PATCH 21/53] [FIX] ruff formatting pt2 --- mne/io/neuralynx/neuralynx.py | 97 ++++++++++++++++++++++------------- 1 file changed, 61 insertions(+), 36 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 9bb07140496..4d267f21e75 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -19,7 +19,6 @@ class AnalogSignalGap(object): """ def __init__(self, signal, units, sampling_rate): - self.signal = signal self.units = units self.sampling_rate = sampling_rate @@ -27,13 +26,14 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], - units=self.units, - sampling_rate=self.sampling_rate) + sig = AnalogSignal( + signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate, + ) return sig - @fill_doc def read_raw_neuralynx( fname, *, preload=False, exclude_fname_patterns=None, verbose=None @@ -127,8 +127,12 @@ def __init__( ) # get segment start/stop times - start_times = np.array([nlx_reader.segment_t_start(block_id, i) for i in range(n_segments)]) - stop_times = np.array([nlx_reader.segment_t_stop(block_id, i) for i in range(n_segments)]) + start_times = np.array( + [nlx_reader.segment_t_start(block_id, i) for i in range(n_segments)] + ) + stop_times = np.array( + [nlx_reader.segment_t_stop(block_id, i) for i in range(n_segments)] + ) # find discontinuous boundaries (of length n-1) next_start_times = start_times[1::] @@ -138,7 +142,7 @@ def __init__( # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) - delta = 1/info["sfreq"] + delta = 1 / info["sfreq"] gaps = seg_diffs > delta has_gaps = gaps.any() @@ -147,16 +151,22 @@ def __init__( gap_annotations = {} if has_gaps: - - logger.info(f"N = {gaps.sum()} discontinuous Neo segments detected with delta > {delta} sec.\n(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})") + logger.info( + f"N = {gaps.sum()} discontinuous Neo segments detected " + + "with delta > {delta} sec. +\n" + + "(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" + ) gap_starts = stop_times[:-1][gaps] # gap starts at segment offset gap_stops = start_times[1::][gaps] # gap stops at segment onset # (n_gaps,) array of ints giving number of samples per inferred gap gap_n_samps = np.array( - [len(np.arange(on_spl, off_spl)) - for on_spl, off_spl in zip(gap_starts*info["sfreq"], gap_stops*info["sfreq"]) + [ + len(np.arange(on_spl, off_spl)) + for on_spl, off_spl in zip( + gap_starts * info["sfreq"], gap_stops * info["sfreq"] + ) ] ) @@ -170,68 +180,82 @@ def __init__( # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( - [np.full(len(start_times), fill_value=0), - np.full(len(gap_starts), fill_value=1) + [ + np.full(len(start_times), fill_value=0), + np.full(len(gap_starts), fill_value=1), ] ) gap_indicator = gap_indicator[all_starts_ids].astype(bool) # store this in a dict to be passed to _raw_extras seg_gap_dict = { - "onsets": all_starts, # onsets in seconds + "onsets": all_starts, # onsets in seconds "offsets": all_stops, "gap_n_samps": gap_n_samps, - "isgap": gap_indicator, # 0 (data segment) or 1 (invalid segment for BAD_SKIP_ACQ) + "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) } # TMP: annotations dict for use with mne.Annotations - gap_annotations = dict(onset=gap_starts, duration=seg_diffs[gaps], orig_time=None, description="BAD_ACQ_SKIP") + gap_annotations = dict( + onset=gap_starts, + duration=seg_diffs[gaps], + orig_time=None, + description="BAD_ACQ_SKIP", + ) - gap_segment_sizes = [ - n for n in gap_n_samps - ] + gap_segment_sizes = [n for n in gap_n_samps] else: - logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") + logger.info( + f"All Neo segments temporally continuous at {delta} sec precision." + ) # check that segment[-1] stop and segment[i] start times # matched to microsecond precision (1e-6) - #breakpoint() - #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, + # breakpoint() + # assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, # err_msg="Segments start/end times are not temporally contiguous." - #) + # ) valid_segment_sizes = [ - nlx_reader.get_signal_size(block_id, i) - for i in range(n_segments) + nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] if has_gaps: - sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[all_starts_ids] + sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[ + all_starts_ids + ] else: sizes_sorted = np.array(valid_segment_sizes) # now construct an (n_samples,) indicator variable sample2segment = np.concatenate( - [np.full(shape=(n,), fill_value=i) - for i, n in enumerate(sizes_sorted)] + [np.full(shape=(n,), fill_value=i) for i, n in enumerate(sizes_sorted)] ) # construct an array of shape (n_total_samples,) indicating # segment membership for each sample - #sample2segment = np.concatenate( + # sample2segment = np.concatenate( # [ # np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) # for i in range(n_segments) # ] - #) + # ) super(RawNeuralynx, self).__init__( info=info, last_samps=[n_total_samples - 1], filenames=[fname], preload=preload, - raw_extras=[dict(smp2seg=sample2segment, exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, gap_annotations=gap_annotations)], + raw_extras=[ + dict( + smp2seg=sample2segment, + exclude_fnames=exclude_fnames, + segment_sizes=sizes_sorted, + seg_gap_dict=seg_gap_dict, + gap_annotations=gap_annotations, + ) + ], ) def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): @@ -253,12 +277,12 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): ) == len(neo_block[0].segments) # collect sizes of each segment - #segment_sizes = np.array( + # segment_sizes = np.array( # [ # nlx_reader.get_signal_size(0, segment_id) # for segment_id in range(len(neo_block[0].segments)) # ] - #) + # ) segment_sizes = self._raw_extras[0]["segment_sizes"] @@ -311,7 +335,6 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): # if gaps were detected correctly insert gap Segments in between valid Segments if self._raw_extras[0]["seg_gap_dict"]: - gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] @@ -320,7 +343,9 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): n_chans = np.arange(idx.start, idx.stop, idx.step).size for seg, n in zip(gap_segments, gap_samples): - asig = AnalogSignalGap(signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz) + asig = AnalogSignalGap( + signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz + ) seg.analogsignals.append(asig) n_total_segments = len(neo_block[0].segments + gap_segments) From b6cfa3ad4b53df6e1f24f4649c1517344c5635b6 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:26:41 -0500 Subject: [PATCH 22/53] [FIX] more style fixes --- mne/io/neuralynx/neuralynx.py | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 4d267f21e75..5022f18be46 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -13,9 +13,11 @@ class AnalogSignalGap(object): - """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """Dummy object to represent gaps in Neuralynx data. + + Creates a AnalogSignalProxy-like object. + Propagate `signal`, `units`, and `sampling_rate` attributes + to the `AnalogSignal` object returned by `load()`. """ def __init__(self, signal, units, sampling_rate): @@ -24,7 +26,7 @@ def __init__(self, signal, units, sampling_rate): self.sampling_rate = sampling_rate def load(self, channel_indexes): - """Dummy method such that it returns object and we access .magnitude""" + """Return AnalogSignal object.""" # self.magnitude = self.magnitude[channel_indexes, :] sig = AnalogSignal( signal=self.signal[channel_indexes, :], @@ -210,13 +212,6 @@ def __init__( f"All Neo segments temporally continuous at {delta} sec precision." ) - # check that segment[-1] stop and segment[i] start times - # matched to microsecond precision (1e-6) - # breakpoint() - # assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, - # err_msg="Segments start/end times are not temporally contiguous." - # ) - valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] From 0f8bf803b984fdc5323090c636156d6952051870 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Thu, 7 Dec 2023 16:27:56 -0500 Subject: [PATCH 23/53] [WIP] prototype to handle gaps --- mne/io/neuralynx/neuralynx.py | 22 ++++++++++++++++++++++ 1 file changed, 22 insertions(+) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 5022f18be46..4a062b8e42b 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -9,9 +9,31 @@ from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose +from ...annotations import Annotations from ..base import BaseRaw +from neo import AnalogSignal + +class AnalogSignalGap(object): + + def __init__(self, signal, units, sampling_rate): + + self.signal = signal + self.units = units + self.sampling_rate = sampling_rate + + def load(self, channel_indexes): + """Dummy method such that it returns object and we access .magnitude""" + + # self.magnitude = self.magnitude[channel_indexes, :] + sig = AnalogSignal(signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate) + return sig + + + class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data. From 11de11002713ab59fc132708df8ab5ec80e93b3c Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:40:27 -0500 Subject: [PATCH 24/53] [ENH] define gap as > than the sampling period --- mne/io/neuralynx/neuralynx.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 4a062b8e42b..686c555f4f4 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -16,7 +16,10 @@ from neo import AnalogSignal class AnalogSignalGap(object): - + """Dummy object to represent gaps in Neuralynx data as + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """ def __init__(self, signal, units, sampling_rate): self.signal = signal From d9dbc1207598160ed34e949e1c8bb7e0fd691961 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 11:42:13 -0500 Subject: [PATCH 25/53] [MAINT] ruff --- mne/io/neuralynx/neuralynx.py | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 686c555f4f4..54b4662de34 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -9,17 +9,15 @@ from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose -from ...annotations import Annotations from ..base import BaseRaw -from neo import AnalogSignal - class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + AnalogSignalProxy-like objects. Propagate `signal`, `units`, and + `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. """ + def __init__(self, signal, units, sampling_rate): self.signal = signal @@ -28,11 +26,10 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" - # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], + sig = AnalogSignal(signal=self.signal[channel_indexes, :], units=self.units, - sampling_rate=self.sampling_rate) + sampling_rate=self.sampling_rate) return sig From ab4f9c0f5b1129a0f2802658d15c265862644e4a Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:19:53 -0500 Subject: [PATCH 26/53] [FIX] ruff formatting pt2 --- mne/io/neuralynx/neuralynx.py | 20 +++++++++++++------- 1 file changed, 13 insertions(+), 7 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 54b4662de34..277c2759836 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -19,7 +19,6 @@ class AnalogSignalGap(object): """ def __init__(self, signal, units, sampling_rate): - self.signal = signal self.units = units self.sampling_rate = sampling_rate @@ -27,9 +26,11 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Dummy method such that it returns object and we access .magnitude""" # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal(signal=self.signal[channel_indexes, :], - units=self.units, - sampling_rate=self.sampling_rate) + sig = AnalogSignal( + signal=self.signal[channel_indexes, :], + units=self.units, + sampling_rate=self.sampling_rate, + ) return sig @@ -230,9 +231,14 @@ def __init__( gap_segment_sizes = [n for n in gap_n_samps] else: - logger.info( - f"All Neo segments temporally continuous at {delta} sec precision." - ) + logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") + + # check that segment[-1] stop and segment[i] start times + # matched to microsecond precision (1e-6) + #breakpoint() + #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, + # err_msg="Segments start/end times are not temporally contiguous." + #) valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) From 2a751c16b706737c1707d5f99bd735f90ba704a4 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Fri, 8 Dec 2023 12:26:41 -0500 Subject: [PATCH 27/53] [FIX] more style fixes --- mne/io/neuralynx/neuralynx.py | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 277c2759836..85943ec0e14 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -13,9 +13,11 @@ class AnalogSignalGap(object): - """Dummy object to represent gaps in Neuralynx data as - AnalogSignalProxy-like objects. Propagate `signal`, `units`, and - `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + """Dummy object to represent gaps in Neuralynx data. + + Creates a AnalogSignalProxy-like object. + Propagate `signal`, `units`, and `sampling_rate` attributes + to the `AnalogSignal` object returned by `load()`. """ def __init__(self, signal, units, sampling_rate): @@ -24,7 +26,7 @@ def __init__(self, signal, units, sampling_rate): self.sampling_rate = sampling_rate def load(self, channel_indexes): - """Dummy method such that it returns object and we access .magnitude""" + """Return AnalogSignal object.""" # self.magnitude = self.magnitude[channel_indexes, :] sig = AnalogSignal( signal=self.signal[channel_indexes, :], @@ -233,13 +235,6 @@ def __init__( else: logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") - # check that segment[-1] stop and segment[i] start times - # matched to microsecond precision (1e-6) - #breakpoint() - #assert_allclose(stop_times[:-1]-start_times[1::], 0, atol=1e-3, - # err_msg="Segments start/end times are not temporally contiguous." - #) - valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] From 9d7c9cd9e87b5ed7ca9863dadbaefb799b6b790d Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 11 Dec 2023 16:21:53 -0500 Subject: [PATCH 28/53] [FIX] add annotations dict to _raw_extras --- mne/io/neuralynx/neuralynx.py | 82 +++++++++++++++++------------------ 1 file changed, 40 insertions(+), 42 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 85943ec0e14..a6af9fe9442 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -18,6 +18,12 @@ class AnalogSignalGap(object): Creates a AnalogSignalProxy-like object. Propagate `signal`, `units`, and `sampling_rate` attributes to the `AnalogSignal` object returned by `load()`. + + Parameters + ---------- + signal : array-like + Array of shape (n_channels, n_samples) containing the data. + """ def __init__(self, signal, units, sampling_rate): @@ -29,7 +35,7 @@ def load(self, channel_indexes): """Return AnalogSignal object.""" # self.magnitude = self.magnitude[channel_indexes, :] sig = AnalogSignal( - signal=self.signal[channel_indexes, :], + signal=self.signal[:, channel_indexes], units=self.units, sampling_rate=self.sampling_rate, ) @@ -148,10 +154,6 @@ def __init__( # the sample sizes of all segments n_segments = nlx_reader.header["nb_segment"][0] block_id = 0 # assumes there's only one block of recording - n_total_samples = sum( - nlx_reader.get_signal_size(block_id, segment) - for segment in range(n_segments) - ) # get segment start/stop times start_times = np.array( @@ -168,20 +170,19 @@ def __init__( # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) - + logger.info("Checking for temporal discontinuities in Neo data segments.") delta = 1 / info["sfreq"] gaps = seg_diffs > delta has_gaps = gaps.any() seg_gap_dict = {} - gap_segment_sizes = [] - gap_annotations = {} if has_gaps: logger.info( - f"N = {gaps.sum()} discontinuous Neo segments detected " + - "with delta > {delta} sec. +\n" + - "(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" + f"N = {gaps.sum()} discontinuous Neo segments detected " + + f"with delta > {delta} sec.\n" + + f"(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" + + "\nAnnotating gaps as BAD_ACQ_SKIP." ) gap_starts = stop_times[:-1][gaps] # gap starts at segment offset @@ -197,13 +198,8 @@ def __init__( ] ) - # add the inferred gaps into the right place in the segment list + # get sort indices for all segments (valid and gap) in ascending order all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) - all_stops_ids = np.argsort(np.concatenate([stop_times, gap_stops])) - - # sort the valid segment and gap times by time - all_starts = np.concatenate([start_times, gap_starts])[all_starts_ids] - all_stops = np.concatenate([stop_times, gap_stops])[all_stops_ids] # variable indicating whether each segment is a gap or not gap_indicator = np.concatenate( @@ -216,22 +212,10 @@ def __init__( # store this in a dict to be passed to _raw_extras seg_gap_dict = { - "onsets": all_starts, # onsets in seconds - "offsets": all_stops, "gap_n_samps": gap_n_samps, "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) } - # TMP: annotations dict for use with mne.Annotations - gap_annotations = dict( - onset=gap_starts, - duration=seg_diffs[gaps], - orig_time=None, - description="BAD_ACQ_SKIP", - ) - - gap_segment_sizes = [n for n in gap_n_samps] - else: logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") @@ -240,7 +224,7 @@ def __init__( ] if has_gaps: - sizes_sorted = np.concatenate([valid_segment_sizes, gap_segment_sizes])[ + sizes_sorted = np.concatenate([valid_segment_sizes, gap_n_samps])[ all_starts_ids ] else: @@ -251,18 +235,30 @@ def __init__( [np.full(shape=(n,), fill_value=i) for i, n in enumerate(sizes_sorted)] ) - # construct an array of shape (n_total_samples,) indicating - # segment membership for each sample - # sample2segment = np.concatenate( - # [ - # np.full(shape=(nlx_reader.get_signal_size(block_id, i),), fill_value=i) - # for i in range(n_segments) - # ] - # ) + if has_gaps: + gap_seg_ids = np.unique(sample2segment)[gap_indicator] + gap_start_ids = np.array( + [np.where(sample2segment == seg_id)[0][0] for seg_id in gap_seg_ids] + ) + + # recreate time axis for gap annotations + mne_times = np.arange(0, len(sample2segment)) / info["sfreq"] + + assert len(gap_start_ids) == len(gap_n_samps) + annotations = dict( + onset=[mne_times[onset_id] for onset_id in gap_start_ids], + duration=[ + mne_times[onset_id + (n - 1)] - mne_times[onset_id] + for onset_id, n in zip(gap_start_ids, gap_n_samps) + ], + descriptions=["BAD_ACQ_SKIP" for _ in gap_start_ids], + ) + else: + annotations = None super(RawNeuralynx, self).__init__( info=info, - last_samps=[n_total_samples - 1], + last_samps=[sizes_sorted.sum() - 1], filenames=[fname], preload=preload, raw_extras=[ @@ -271,7 +267,7 @@ def __init__( exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, - gap_annotations=gap_annotations, + gap_annotations=annotations, ) ], ) @@ -356,7 +352,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] - # create AnalogSignal objects representing gap data + # create AnalogSignal objects representing gap data filled with 0's sfreq = nlx_reader.get_signal_sampling_rate() n_chans = np.arange(idx.start, idx.stop, idx.step).size @@ -368,12 +364,14 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): n_total_segments = len(neo_block[0].segments + gap_segments) segments_arr = np.zeros((n_total_segments,), dtype=object) + + # insert inferred gap segments at the right place in between valid segments isgap = self._raw_extras[0]["seg_gap_dict"]["isgap"] segments_arr[~isgap] = neo_block[0].segments segments_arr[isgap] = gap_segments # now load data from selected segments/channels via - # neo.Segment.AnalogSignal.load() + # neo.Segment.AnalogSignal.load() or AnalogSignalGap.load() all_data = np.concatenate( [ signal.load(channel_indexes=idx).magnitude[ From 395b1ab777437abb45556d744625cd689d983c68 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 11 Dec 2023 16:46:26 -0500 Subject: [PATCH 29/53] [MAINT] remove commented section --- mne/io/neuralynx/neuralynx.py | 8 -------- 1 file changed, 8 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index a6af9fe9442..ab24ecc2dd7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -290,14 +290,6 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): [len(segment.analogsignals) for segment in neo_block[0].segments] ) == len(neo_block[0].segments) - # collect sizes of each segment - # segment_sizes = np.array( - # [ - # nlx_reader.get_signal_size(0, segment_id) - # for segment_id in range(len(neo_block[0].segments)) - # ] - # ) - segment_sizes = self._raw_extras[0]["segment_sizes"] # construct a (n_segments, 2) array of the first and last From 98cef5aa62c8db0cb9e96ef563d45b3452d995d0 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 08:58:24 -0500 Subject: [PATCH 30/53] [FIX] remove doubled code from conflicts --- mne/io/neuralynx/neuralynx.py | 29 +++-------------------------- 1 file changed, 3 insertions(+), 26 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index ab24ecc2dd7..476a319a6b7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -42,31 +42,6 @@ def load(self, channel_indexes): return sig - -class AnalogSignalGap(object): - """Dummy object to represent gaps in Neuralynx data. - - Creates a AnalogSignalProxy-like object. - Propagate `signal`, `units`, and `sampling_rate` attributes - to the `AnalogSignal` object returned by `load()`. - """ - - def __init__(self, signal, units, sampling_rate): - self.signal = signal - self.units = units - self.sampling_rate = sampling_rate - - def load(self, channel_indexes): - """Return AnalogSignal object.""" - # self.magnitude = self.magnitude[channel_indexes, :] - sig = AnalogSignal( - signal=self.signal[channel_indexes, :], - units=self.units, - sampling_rate=self.sampling_rate, - ) - return sig - - @fill_doc def read_raw_neuralynx( fname, *, preload=False, exclude_fname_patterns=None, verbose=None @@ -217,7 +192,9 @@ def __init__( } else: - logger.info(f"All Neo segments temporally continuous at {delta} sec precision.") + logger.info( + f"All Neo segments temporally continuous at {delta} sec precision." + ) valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) From 7902bd806c1d57de692fac0ace44ab778d5f962e Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 09:02:12 -0500 Subject: [PATCH 31/53] [MAINT] fill docstring --- mne/io/neuralynx/neuralynx.py | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 476a319a6b7..81cc3f6f80a 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -17,13 +17,21 @@ class AnalogSignalGap(object): Creates a AnalogSignalProxy-like object. Propagate `signal`, `units`, and `sampling_rate` attributes - to the `AnalogSignal` object returned by `load()`. + to the `AnalogSignal` init returned by `load()`. Parameters ---------- signal : array-like Array of shape (n_channels, n_samples) containing the data. + units : str + Units of the data. (e.g., 'uV') + sampling_rate : quantity + Sampling rate of the data. (e.g., 4000 * pq.Hz) + Returns + ------- + sig : instance of AnalogSignal + A AnalogSignal object representing a gap in Neuralynx data. """ def __init__(self, signal, units, sampling_rate): From fce464720809a9f8b5cba8598fc366b36d1028b0 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 09:45:58 -0500 Subject: [PATCH 32/53] [FIX] move neo imports at the top --- mne/io/neuralynx/neuralynx.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 81cc3f6f80a..c825f8195f5 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -4,13 +4,16 @@ import os import numpy as np -from neo import AnalogSignal from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose from ..base import BaseRaw +_soft_import("neo", "Reading NeuralynxIO files", strict=True) +from neo import AnalogSignal +from neo.io import NeuralynxIO + class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data. @@ -98,9 +101,6 @@ def __init__( exclude_fname_patterns=None, verbose=None, ): - _soft_import("neo", "Reading NeuralynxIO files", strict=True) - from neo.io import NeuralynxIO - fname = _check_fname(fname, "read", True, "fname", need_dir=True) logger.info(f"Checking files in {fname}") From e895c0b8a63aa950af254282c610a136e3c6418d Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 09:51:36 -0500 Subject: [PATCH 33/53] [FIX] move neo soft imports within methods --- mne/io/neuralynx/neuralynx.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index c825f8195f5..2d51a13ed08 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -10,10 +10,6 @@ from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose from ..base import BaseRaw -_soft_import("neo", "Reading NeuralynxIO files", strict=True) -from neo import AnalogSignal -from neo.io import NeuralynxIO - class AnalogSignalGap(object): """Dummy object to represent gaps in Neuralynx data. @@ -44,7 +40,9 @@ def __init__(self, signal, units, sampling_rate): def load(self, channel_indexes): """Return AnalogSignal object.""" - # self.magnitude = self.magnitude[channel_indexes, :] + _soft_import("neo", "Reading NeuralynxIO files", strict=True) + from neo import AnalogSignal + sig = AnalogSignal( signal=self.signal[:, channel_indexes], units=self.units, @@ -103,6 +101,9 @@ def __init__( ): fname = _check_fname(fname, "read", True, "fname", need_dir=True) + _soft_import("neo", "Reading NeuralynxIO files", strict=True) + from neo.io import NeuralynxIO + logger.info(f"Checking files in {fname}") # construct a list of filenames to ignore From 7f5bfe736b2ce78a2435d6e88e6040c13143e674 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Tue, 12 Dec 2023 16:30:22 -0500 Subject: [PATCH 34/53] [ENH] set annotations post init --- mne/io/neuralynx/neuralynx.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 2d51a13ed08..7da2ccda089 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -7,6 +7,7 @@ from ..._fiff.meas_info import create_info from ..._fiff.utils import _mult_cal_one +from ...annotations import Annotations from ...utils import _check_fname, _soft_import, fill_doc, logger, verbose from ..base import BaseRaw @@ -221,6 +222,7 @@ def __init__( [np.full(shape=(n,), fill_value=i) for i, n in enumerate(sizes_sorted)] ) + # construct Annotations() if has_gaps: gap_seg_ids = np.unique(sample2segment)[gap_indicator] gap_start_ids = np.array( @@ -231,13 +233,13 @@ def __init__( mne_times = np.arange(0, len(sample2segment)) / info["sfreq"] assert len(gap_start_ids) == len(gap_n_samps) - annotations = dict( + annotations = Annotations( onset=[mne_times[onset_id] for onset_id in gap_start_ids], duration=[ mne_times[onset_id + (n - 1)] - mne_times[onset_id] for onset_id, n in zip(gap_start_ids, gap_n_samps) ], - descriptions=["BAD_ACQ_SKIP" for _ in gap_start_ids], + description=["BAD_ACQ_SKIP" for _ in gap_start_ids], ) else: annotations = None @@ -253,11 +255,12 @@ def __init__( exclude_fnames=exclude_fnames, segment_sizes=sizes_sorted, seg_gap_dict=seg_gap_dict, - gap_annotations=annotations, ) ], ) + self.set_annotations(annotations) + def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" from neo import Segment From 730cbf7e6602e2863e9e613542ede6493b5d7a25 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Wed, 13 Dec 2023 15:31:52 -0500 Subject: [PATCH 35/53] [MAINT] comment for quantities dependency Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 1 + 1 file changed, 1 insertion(+) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 7da2ccda089..607a6b5bb82 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -265,6 +265,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" from neo import Segment from neo.io import NeuralynxIO + # quantities is a dependency of neo so we are guaranteed it exists from quantities import Hz nlx_reader = NeuralynxIO( From 37197c7b18546bb77af7e075d691594d8ee1e854 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Wed, 13 Dec 2023 17:29:10 -0500 Subject: [PATCH 36/53] [ENH] remove branching code path --- mne/io/neuralynx/neuralynx.py | 160 ++++++++++++++++------------------ 1 file changed, 77 insertions(+), 83 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 607a6b5bb82..04a6a5ee243 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -156,66 +156,57 @@ def __init__( # mark as discontinuous any two segments that have # start/stop delta larger than sampling period (1/sampling_rate) logger.info("Checking for temporal discontinuities in Neo data segments.") - delta = 1 / info["sfreq"] + delta = 1.5 / info["sfreq"] gaps = seg_diffs > delta - has_gaps = gaps.any() seg_gap_dict = {} - if has_gaps: - logger.info( - f"N = {gaps.sum()} discontinuous Neo segments detected " - + f"with delta > {delta} sec.\n" - + f"(max = {seg_diffs[gaps].max()} sec, min = {seg_diffs[gaps].min()})" - + "\nAnnotating gaps as BAD_ACQ_SKIP." - ) + logger.info( + f"N = {gaps.sum()} discontinuous Neo segments detected " + + f"with delta > {delta} sec." + + "Annotating gaps as BAD_ACQ_SKIP." + if gaps.any() + else "No discontinuities detected." + ) - gap_starts = stop_times[:-1][gaps] # gap starts at segment offset - gap_stops = start_times[1::][gaps] # gap stops at segment onset + gap_starts = stop_times[:-1][gaps] # gap starts at segment offset + gap_stops = start_times[1::][gaps] # gap stops at segment onset - # (n_gaps,) array of ints giving number of samples per inferred gap - gap_n_samps = np.array( - [ - len(np.arange(on_spl, off_spl)) - for on_spl, off_spl in zip( - gap_starts * info["sfreq"], gap_stops * info["sfreq"] - ) - ] - ) - - # get sort indices for all segments (valid and gap) in ascending order - all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) + # (n_gaps,) array of ints giving number of samples per inferred gap + gap_n_samps = np.array( + [ + len(np.arange(on_spl, off_spl)) + for on_spl, off_spl in zip( + gap_starts * info["sfreq"], gap_stops * info["sfreq"] + ) + ] + ).astype(int) # force an int array (if no gaps, empty array is a float) - # variable indicating whether each segment is a gap or not - gap_indicator = np.concatenate( - [ - np.full(len(start_times), fill_value=0), - np.full(len(gap_starts), fill_value=1), - ] - ) - gap_indicator = gap_indicator[all_starts_ids].astype(bool) + # get sort indices for all segments (valid and gap) in ascending order + all_starts_ids = np.argsort(np.concatenate([start_times, gap_starts])) - # store this in a dict to be passed to _raw_extras - seg_gap_dict = { - "gap_n_samps": gap_n_samps, - "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) - } + # variable indicating whether each segment is a gap or not + gap_indicator = np.concatenate( + [ + np.full(len(start_times), fill_value=0), + np.full(len(gap_starts), fill_value=1), + ] + ) + gap_indicator = gap_indicator[all_starts_ids].astype(bool) - else: - logger.info( - f"All Neo segments temporally continuous at {delta} sec precision." - ) + # store this in a dict to be passed to _raw_extras + seg_gap_dict = { + "gap_n_samps": gap_n_samps, + "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) + } valid_segment_sizes = [ nlx_reader.get_signal_size(block_id, i) for i in range(n_segments) ] - if has_gaps: - sizes_sorted = np.concatenate([valid_segment_sizes, gap_n_samps])[ - all_starts_ids - ] - else: - sizes_sorted = np.array(valid_segment_sizes) + sizes_sorted = np.concatenate([valid_segment_sizes, gap_n_samps])[ + all_starts_ids + ] # now construct an (n_samples,) indicator variable sample2segment = np.concatenate( @@ -223,26 +214,23 @@ def __init__( ) # construct Annotations() - if has_gaps: - gap_seg_ids = np.unique(sample2segment)[gap_indicator] - gap_start_ids = np.array( - [np.where(sample2segment == seg_id)[0][0] for seg_id in gap_seg_ids] - ) + gap_seg_ids = np.unique(sample2segment)[gap_indicator] + gap_start_ids = np.array( + [np.where(sample2segment == seg_id)[0][0] for seg_id in gap_seg_ids] + ) - # recreate time axis for gap annotations - mne_times = np.arange(0, len(sample2segment)) / info["sfreq"] - - assert len(gap_start_ids) == len(gap_n_samps) - annotations = Annotations( - onset=[mne_times[onset_id] for onset_id in gap_start_ids], - duration=[ - mne_times[onset_id + (n - 1)] - mne_times[onset_id] - for onset_id, n in zip(gap_start_ids, gap_n_samps) - ], - description=["BAD_ACQ_SKIP" for _ in gap_start_ids], - ) - else: - annotations = None + # recreate time axis for gap annotations + mne_times = np.arange(0, len(sample2segment)) / info["sfreq"] + + assert len(gap_start_ids) == len(gap_n_samps) + annotations = Annotations( + onset=[mne_times[onset_id] for onset_id in gap_start_ids], + duration=[ + mne_times[onset_id + (n - 1)] - mne_times[onset_id] + for onset_id, n in zip(gap_start_ids, gap_n_samps) + ], + description=["BAD_ACQ_SKIP" for _ in gap_start_ids], + ) super(RawNeuralynx, self).__init__( info=info, @@ -265,6 +253,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): """Read a chunk of raw data.""" from neo import Segment from neo.io import NeuralynxIO + # quantities is a dependency of neo so we are guaranteed it exists from quantities import Hz @@ -330,27 +319,32 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): segments_arr = np.array(neo_block[0].segments, dtype=object) # if gaps were detected correctly insert gap Segments in between valid Segments - if self._raw_extras[0]["seg_gap_dict"]: - gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] - gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] - - # create AnalogSignal objects representing gap data filled with 0's - sfreq = nlx_reader.get_signal_sampling_rate() - n_chans = np.arange(idx.start, idx.stop, idx.step).size + # if self._raw_extras[0]["seg_gap_dict"]: + gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] + gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] + + # create AnalogSignal objects representing gap data filled with 0's + # breakpoint() + sfreq = nlx_reader.get_signal_sampling_rate() + n_chans = ( + np.arange(idx.start, idx.stop, idx.step).size + if type(idx) is slice + else len(idx) # idx can be a slice or an np.array so check both + ) - for seg, n in zip(gap_segments, gap_samples): - asig = AnalogSignalGap( - signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz - ) - seg.analogsignals.append(asig) + for seg, n in zip(gap_segments, gap_samples): + asig = AnalogSignalGap( + signal=np.zeros((n, n_chans)), units="uV", sampling_rate=sfreq * Hz + ) + seg.analogsignals.append(asig) - n_total_segments = len(neo_block[0].segments + gap_segments) - segments_arr = np.zeros((n_total_segments,), dtype=object) + n_total_segments = len(neo_block[0].segments + gap_segments) + segments_arr = np.zeros((n_total_segments,), dtype=object) - # insert inferred gap segments at the right place in between valid segments - isgap = self._raw_extras[0]["seg_gap_dict"]["isgap"] - segments_arr[~isgap] = neo_block[0].segments - segments_arr[isgap] = gap_segments + # insert inferred gap segments at the right place in between valid segments + isgap = self._raw_extras[0]["seg_gap_dict"]["isgap"] + segments_arr[~isgap] = neo_block[0].segments + segments_arr[isgap] = gap_segments # now load data from selected segments/channels via # neo.Segment.AnalogSignal.load() or AnalogSignalGap.load() From 2c0861c03ecd79b1e65115440318929c4f3acdb3 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Wed, 13 Dec 2023 17:30:57 -0500 Subject: [PATCH 37/53] [ENH] remove np.arange to infer n samples --- mne/io/neuralynx/neuralynx.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 04a6a5ee243..e7da23bc0e7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -175,10 +175,8 @@ def __init__( # (n_gaps,) array of ints giving number of samples per inferred gap gap_n_samps = np.array( [ - len(np.arange(on_spl, off_spl)) - for on_spl, off_spl in zip( - gap_starts * info["sfreq"], gap_stops * info["sfreq"] - ) + int(round(stop * info["sfreq"])) - int(round(start * info["sfreq"])) + for start, stop in zip(gap_starts, gap_stops) ] ).astype(int) # force an int array (if no gaps, empty array is a float) From ab00f83b24c91de8b3f85dae412ed3179e5e4453 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Thu, 14 Dec 2023 09:10:04 -0500 Subject: [PATCH 38/53] [ENH] improve error handling, comment cleanup --- mne/io/neuralynx/neuralynx.py | 22 +++++++++++++--------- 1 file changed, 13 insertions(+), 9 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index e7da23bc0e7..26ca9ed8ad7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -122,12 +122,18 @@ def __init__( try: nlx_reader = NeuralynxIO(dirname=fname, exclude_filename=exclude_fnames) except ValueError as e: - raise ValueError( - "It seems some .ncs channels might have different number of samples. " - + "This is likely due to different sampling rates. " - + "Try excluding them with `exclude_fname_patterns` input arg." - + f"\nOriginal neo.NeuralynxIO.parse_header() ValueError:\n{e}" - ) + # give a more informative error message and what the user can do about it + if "Incompatible section structures across streams" in e.__str__(): + raise ValueError( + "It seems .ncs channels have different numbers of samples. " + + "This is likely due to different sampling rates. " + + "Try reading in only channels with uniform sampling rate " + + "by excluding other channels with `exclude_fname_patterns` " + + "input argument." + + f"\nOriginal neo.NeuralynxRawIO._parse_header() ValueError:\n{e}" + ) + else: + raise e # any other Neo error, just raise it info = create_info( ch_types="seeg", @@ -316,13 +322,11 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): # array containing Segments segments_arr = np.array(neo_block[0].segments, dtype=object) - # if gaps were detected correctly insert gap Segments in between valid Segments - # if self._raw_extras[0]["seg_gap_dict"]: + # if gaps were detected, correctly insert gap Segments in between valid Segments gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] # create AnalogSignal objects representing gap data filled with 0's - # breakpoint() sfreq = nlx_reader.get_signal_sampling_rate() n_chans = ( np.arange(idx.start, idx.stop, idx.step).size From ffb3e1acb8a65d49a6d37d1404ec78c2ba7f44ae Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 11:48:09 -0500 Subject: [PATCH 39/53] [ENH] use str() instea of __str__() Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index d5d67d68130..b3ed72eddfb 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -123,7 +123,7 @@ def __init__( nlx_reader = NeuralynxIO(dirname=fname, exclude_filename=exclude_fnames) except ValueError as e: # give a more informative error message and what the user can do about it - if "Incompatible section structures across streams" in e.__str__(): + if "Incompatible section structures across streams" in str(e): raise ValueError( "It seems .ncs channels have different numbers of samples. " + "This is likely due to different sampling rates. " From 7b7231ba001865e7b32cb29e83b1f387694f484f Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 11:50:22 -0500 Subject: [PATCH 40/53] [ENH] cleaner error traceback Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index b3ed72eddfb..5c36cadb97d 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -133,7 +133,7 @@ def __init__( + f"\nOriginal neo.NeuralynxRawIO._parse_header() ValueError:\n{e}" ) else: - raise e # any other Neo error, just raise it + raise info = create_info( ch_types="seeg", From 2f821da7f14b8a9a5158a9dd264367d107216187 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 11:51:44 -0500 Subject: [PATCH 41/53] [ENH] avoid duplicate info in error traceback Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 5c36cadb97d..25332c185e7 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -131,7 +131,7 @@ def __init__( + "by excluding other channels with `exclude_fname_patterns` " + "input argument." + f"\nOriginal neo.NeuralynxRawIO._parse_header() ValueError:\n{e}" - ) + ) from None else: raise From 2a971bf349a7cb08e3e5ed1a81568b68f236edac Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 11:52:31 -0500 Subject: [PATCH 42/53] [ENH] don't ref private methods Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 25332c185e7..e6833f42a6c 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -130,7 +130,7 @@ def __init__( + "Try reading in only channels with uniform sampling rate " + "by excluding other channels with `exclude_fname_patterns` " + "input argument." - + f"\nOriginal neo.NeuralynxRawIO._parse_header() ValueError:\n{e}" + + f"\nOriginal neo.NeuralynxRawIO ValueError:\n{e}" ) from None else: raise From e496e1269f6c15d1bc82aeb6756010fc2cf3c2f8 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 11:57:02 -0500 Subject: [PATCH 43/53] Apply suggestions from code review Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index e6833f42a6c..0329ae7ecfb 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -201,7 +201,7 @@ def __init__( # store this in a dict to be passed to _raw_extras seg_gap_dict = { "gap_n_samps": gap_n_samps, - "isgap": gap_indicator, # 0 (data segment) or 1 (gap segment) + "isgap": gap_indicator, # False (data segment) or True (gap segment) } valid_segment_sizes = [ @@ -233,7 +233,7 @@ def __init__( mne_times[onset_id + (n - 1)] - mne_times[onset_id] for onset_id, n in zip(gap_start_ids, gap_n_samps) ], - description=["BAD_ACQ_SKIP" for _ in gap_start_ids], + description=["BAD_ACQ_SKIP"] * len(gap_start_ids), ) super(RawNeuralynx, self).__init__( From 2fe302d17658e2a2d5a99fa6efc7033b925d6a38 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 11:58:02 -0500 Subject: [PATCH 44/53] [ENH] index raw_extras with fi Co-authored-by: Eric Larson --- mne/io/neuralynx/neuralynx.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 0329ae7ecfb..21457019515 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -273,7 +273,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): [len(segment.analogsignals) for segment in neo_block[0].segments] ) == len(neo_block[0].segments) - segment_sizes = self._raw_extras[0]["segment_sizes"] + segment_sizes = self._raw_extras[fi]["segment_sizes"] # construct a (n_segments, 2) array of the first and last # sample index for each segment relative to the start of the recording @@ -323,7 +323,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult): segments_arr = np.array(neo_block[0].segments, dtype=object) # if gaps were detected, correctly insert gap Segments in between valid Segments - gap_samples = self._raw_extras[0]["seg_gap_dict"]["gap_n_samps"] + gap_samples = self._raw_extras[fi]["seg_gap_dict"]["gap_n_samps"] gap_segments = [Segment(f"gap-{i}") for i in range(len(gap_samples))] # create AnalogSignal objects representing gap data filled with 0's From eccb4817acd0d4f0b3ba63f977be5a228144329e Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 15:20:27 -0500 Subject: [PATCH 45/53] [ENH] update version and hash --- mne/datasets/config.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/mne/datasets/config.py b/mne/datasets/config.py index b548f5273f2..b7780778f24 100644 --- a/mne/datasets/config.py +++ b/mne/datasets/config.py @@ -88,7 +88,7 @@ # respective repos, and make a new release of the dataset on GitHub. Then # update the checksum in the MNE_DATASETS dict below, and change version # here: ↓↓↓↓↓ ↓↓↓ -RELEASES = dict(testing="0.150", misc="0.27") +RELEASES = dict(testing="0.151", misc="0.27") TESTING_VERSIONED = f'mne-testing-data-{RELEASES["testing"]}' MISC_VERSIONED = f'mne-misc-data-{RELEASES["misc"]}' @@ -112,7 +112,7 @@ # Testing and misc are at the top as they're updated most often MNE_DATASETS["testing"] = dict( archive_name=f"{TESTING_VERSIONED}.tar.gz", - hash="md5:0b7452daef4d19132505b5639d695628", + hash="md5:5832b4d44f0423d22305fa61cb75bc25", url=( "https://codeload.github.com/mne-tools/mne-testing-data/" f'tar.gz/{RELEASES["testing"]}' From cc56bcb8afa08b172078269215dcc1ae5f40df12 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 15:21:25 -0500 Subject: [PATCH 46/53] [MAINT] blank space in info msg --- mne/io/neuralynx/neuralynx.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne/io/neuralynx/neuralynx.py b/mne/io/neuralynx/neuralynx.py index 21457019515..1c007ba5787 100644 --- a/mne/io/neuralynx/neuralynx.py +++ b/mne/io/neuralynx/neuralynx.py @@ -169,7 +169,7 @@ def __init__( logger.info( f"N = {gaps.sum()} discontinuous Neo segments detected " - + f"with delta > {delta} sec." + + f"with delta > {delta} sec. " + "Annotating gaps as BAD_ACQ_SKIP." if gaps.any() else "No discontinuities detected." From e3d4d8ebf0e862a2c3173fb40b79c93a1ece595e Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 15:22:52 -0500 Subject: [PATCH 47/53] [ENH] add test_neuralynx_gaps() --- mne/io/neuralynx/tests/test_neuralynx.py | 105 ++++++++++++++++++++--- 1 file changed, 93 insertions(+), 12 deletions(-) diff --git a/mne/io/neuralynx/tests/test_neuralynx.py b/mne/io/neuralynx/tests/test_neuralynx.py index 21cb73927a8..26b372efc2f 100644 --- a/mne/io/neuralynx/tests/test_neuralynx.py +++ b/mne/io/neuralynx/tests/test_neuralynx.py @@ -15,7 +15,6 @@ testing_path = data_path(download=False) / "neuralynx" - def _nlxheader_to_dict(matdict: Dict) -> Dict: """Convert the read-in "Header" field into a dict. @@ -65,14 +64,42 @@ def _read_nlx_mat_chan(matfile: str) -> np.ndarray: return x -mne_testing_ncs = [ - "LAHC1.ncs", - "LAHC2.ncs", - "LAHC3.ncs", - "LAHCu1.ncs", # the 'u' files are going to be filtered out - "xAIR1.ncs", - "xEKG1.ncs", -] +def _read_nlx_mat_chan_keep_gaps(matfile: str) -> np.ndarray: + """Read a single channel from a Neuralynx .mat file.""" + mat = loadmat(matfile) + + hdr_dict = _nlxheader_to_dict(mat) + + # Nlx2MatCSC.m reads the data in N equal-sized (512-item) chunks + # this array (1, n_chunks) stores the number of valid samples + # per chunk (the last chunk is usually shorter) + n_valid_samples = mat["NumberOfValidSamples"].ravel() + + # read in the artificial zeros so that + # we can compare with the mne padded arrays + ncs_records_with_gaps = [9, 15, 20] + for i in ncs_records_with_gaps: + n_valid_samples[i] = 512 + + # concatenate chunks, respecting the number of valid samples + x = np.concatenate( + [mat["Samples"][0:n, i] for i, n in enumerate(n_valid_samples)] + ) # in ADBits + + # this value is the same for all channels and + # converts data from ADBits to Volts + conversionf = literal_eval(hdr_dict["ADBitVolts"]) + x = x * conversionf + + # if header says input was inverted at acquisition + # (possibly for spike detection or so?), flip it back + # NeuralynxIO does this under the hood in NeuralynxIO.parse_header() + # see this discussion: https://github.com/NeuralEnsemble/python-neo/issues/819 + if hdr_dict["InputInverted"] == "True": + x *= -1 + + return x + expected_chan_names = ["LAHC1", "LAHC2", "LAHC3", "xAIR1", "xEKG1"] @@ -84,11 +111,18 @@ def test_neuralynx(): from neo.io import NeuralynxIO - excluded_ncs_files = ["LAHCu1.ncs", "LAHCu2.ncs", "LAHCu3.ncs"] + excluded_ncs_files = [ + "LAHCu1.ncs", + "LAHC1_3_gaps.ncs", + "LAHC2_3_gaps.ncs", + ] # ==== MNE-Python ==== # + fname_patterns = ["*u*.ncs", "*3_gaps.ncs"] raw = read_raw_neuralynx( - fname=testing_path, preload=True, exclude_fname_patterns=["*u*.ncs"] + fname=testing_path, + preload=True, + exclude_fname_patterns=fname_patterns, ) # test that channel selection worked @@ -136,5 +170,52 @@ def test_neuralynx(): ) # data _test_raw_reader( - read_raw_neuralynx, fname=testing_path, exclude_fname_patterns=["*u*.ncs"] + read_raw_neuralynx, + fname=testing_path, + exclude_fname_patterns=fname_patterns, + ) + + +@requires_testing_data +def test_neuralynx_gaps(): + """Test gap detection.""" + # ignore files with no gaps + ignored_ncs_files = [ + "LAHC1.ncs", + "LAHC2.ncs", + "LAHC3.ncs", + "xAIR1.ncs", + "xEKG1.ncs", + "LAHCu1.ncs", + ] + raw = read_raw_neuralynx( + fname=testing_path, + preload=True, + exclude_fname_patterns=ignored_ncs_files, + ) + mne_y, _ = raw.get_data(return_times=True) # in V + + # there should be 2 channels with 3 gaps (of 130 samples in total) + n_expected_gaps = 3 + n_expected_missing_samples = 130 + assert len(raw.annotations) == n_expected_gaps, "Wrong number of gaps detected" + assert ( + (mne_y[0, :] == 0).sum() == n_expected_missing_samples + ), "Number of true and inferred missing samples differ" + + # read in .mat files containing original gaps + matchans = ["LAHC1_3_gaps.mat", "LAHC2_3_gaps.mat"] + + # (n_chan, n_samples) array, in V + mat_y = np.stack( + [ + _read_nlx_mat_chan_keep_gaps(os.path.join(testing_path, ch)) + for ch in matchans + ] + ) + + # compare originally modified .ncs arrays with MNE-padded arrays + # and test that we back-inserted 0's at the right places + assert_allclose( + mne_y, mat_y, rtol=1e-6, err_msg="MNE and Nlx2MatCSC.m not all close" ) From 81ec4da9b35dfc3e679cc340473793592ba8ed3c Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 15:29:41 -0500 Subject: [PATCH 48/53] [FIX] formater --- mne/io/neuralynx/tests/test_neuralynx.py | 1 + 1 file changed, 1 insertion(+) diff --git a/mne/io/neuralynx/tests/test_neuralynx.py b/mne/io/neuralynx/tests/test_neuralynx.py index 26b372efc2f..65c5b7f5d12 100644 --- a/mne/io/neuralynx/tests/test_neuralynx.py +++ b/mne/io/neuralynx/tests/test_neuralynx.py @@ -15,6 +15,7 @@ testing_path = data_path(download=False) / "neuralynx" + def _nlxheader_to_dict(matdict: Dict) -> Dict: """Convert the read-in "Header" field into a dict. From c26c7ea45d6a523956e2becb0fd1aee186bc7250 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 15:31:28 -0500 Subject: [PATCH 49/53] [MAINT] docstring update --- mne/io/neuralynx/tests/test_neuralynx.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne/io/neuralynx/tests/test_neuralynx.py b/mne/io/neuralynx/tests/test_neuralynx.py index 65c5b7f5d12..53b71886da4 100644 --- a/mne/io/neuralynx/tests/test_neuralynx.py +++ b/mne/io/neuralynx/tests/test_neuralynx.py @@ -66,7 +66,7 @@ def _read_nlx_mat_chan(matfile: str) -> np.ndarray: def _read_nlx_mat_chan_keep_gaps(matfile: str) -> np.ndarray: - """Read a single channel from a Neuralynx .mat file.""" + """Read a single channel from a Neuralynx .mat file and keep invalid samples.""" mat = loadmat(matfile) hdr_dict = _nlxheader_to_dict(mat) From dccc6528d1a72fc195cdd3e293b2fb0be6bd4fc5 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 15:32:18 -0500 Subject: [PATCH 50/53] [MAINT] update changelog in devel.rst --- doc/changes/devel.rst | 1 + 1 file changed, 1 insertion(+) diff --git a/doc/changes/devel.rst b/doc/changes/devel.rst index d993f4cc26c..cae3e95a065 100644 --- a/doc/changes/devel.rst +++ b/doc/changes/devel.rst @@ -49,6 +49,7 @@ Bugs - Remove incorrect type hints in :func:`mne.io.read_raw_neuralynx` (:gh:`12236` by `Richard Höchenberger`_) - Fix bug where parent directory existence was not checked properly in :meth:`mne.io.Raw.save` (:gh:`12282` by `Eric Larson`_) - ``defusedxml`` is now an optional (rather than required) dependency and needed when reading EGI-MFF data, NEDF data, and BrainVision montages (:gh:`12264` by `Eric Larson`_) +- Correctly handle temporal gaps in Neuralynx `.ncs` files via :func:`mne.io.read_raw_neuralynx` (:gh:`12279` by `Kristijan Armeni`_ and `Eric Larson`_) API changes ~~~~~~~~~~~ From c081277fe2dc33ccba2d323d3e90f79be66b72f0 Mon Sep 17 00:00:00 2001 From: Kristijan Armeni Date: Mon, 18 Dec 2023 16:01:11 -0500 Subject: [PATCH 51/53] [FIX] formatting in changelog text --- doc/changes/devel.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/changes/devel.rst b/doc/changes/devel.rst index cae3e95a065..858ccf019fd 100644 --- a/doc/changes/devel.rst +++ b/doc/changes/devel.rst @@ -49,7 +49,7 @@ Bugs - Remove incorrect type hints in :func:`mne.io.read_raw_neuralynx` (:gh:`12236` by `Richard Höchenberger`_) - Fix bug where parent directory existence was not checked properly in :meth:`mne.io.Raw.save` (:gh:`12282` by `Eric Larson`_) - ``defusedxml`` is now an optional (rather than required) dependency and needed when reading EGI-MFF data, NEDF data, and BrainVision montages (:gh:`12264` by `Eric Larson`_) -- Correctly handle temporal gaps in Neuralynx `.ncs` files via :func:`mne.io.read_raw_neuralynx` (:gh:`12279` by `Kristijan Armeni`_ and `Eric Larson`_) +- Correctly handle temporal gaps in Neuralynx .ncs files via :func:`mne.io.read_raw_neuralynx` (:gh:`12279` by `Kristijan Armeni`_ and `Eric Larson`_) API changes ~~~~~~~~~~~ From b26dae83b016d7159b3e85118b27fa0a26d1a1a9 Mon Sep 17 00:00:00 2001 From: Eric Larson Date: Tue, 19 Dec 2023 11:09:03 -0500 Subject: [PATCH 52/53] FIX: Module level skip and install --- environment.yml | 1 + mne/io/neuralynx/tests/test_neuralynx.py | 4 ++-- mne/utils/config.py | 1 + pyproject.toml | 2 ++ 4 files changed, 6 insertions(+), 2 deletions(-) diff --git a/environment.yml b/environment.yml index 8978dfc64e8..abdb2740444 100644 --- a/environment.yml +++ b/environment.yml @@ -61,3 +61,4 @@ dependencies: - mamba - lazy_loader - defusedxml + - neo diff --git a/mne/io/neuralynx/tests/test_neuralynx.py b/mne/io/neuralynx/tests/test_neuralynx.py index 53b71886da4..1532845ab7a 100644 --- a/mne/io/neuralynx/tests/test_neuralynx.py +++ b/mne/io/neuralynx/tests/test_neuralynx.py @@ -15,6 +15,8 @@ testing_path = data_path(download=False) / "neuralynx" +pytest.importorskip("neo") + def _nlxheader_to_dict(matdict: Dict) -> Dict: """Convert the read-in "Header" field into a dict. @@ -108,8 +110,6 @@ def _read_nlx_mat_chan_keep_gaps(matfile: str) -> np.ndarray: @requires_testing_data def test_neuralynx(): """Test basic reading.""" - pytest.importorskip("neo") - from neo.io import NeuralynxIO excluded_ncs_files = [ diff --git a/mne/utils/config.py b/mne/utils/config.py index 77b94508114..62b4d053012 100644 --- a/mne/utils/config.py +++ b/mne/utils/config.py @@ -684,6 +684,7 @@ def sys_info( "mne-connectivity", "mne-icalabel", "mne-bids-pipeline", + "neo", "", ) if dependencies == "developer": diff --git a/pyproject.toml b/pyproject.toml index 39c6876e43d..092e4dd102a 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -104,6 +104,7 @@ full = [ "pybv", "snirf", "defusedxml", + "neo", ] # Dependencies for running the test infrastructure @@ -135,6 +136,7 @@ test_extra = [ "imageio>=2.6.1", "imageio-ffmpeg>=0.4.1", "snirf", + "neo", ] # Dependencies for building the docuemntation From 965a30d9749dc6f0a7f83b2e4cbe686112d1116d Mon Sep 17 00:00:00 2001 From: Eric Larson Date: Tue, 19 Dec 2023 11:32:17 -0500 Subject: [PATCH 53/53] FIX: Name --- environment.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/environment.yml b/environment.yml index abdb2740444..96c89fe472b 100644 --- a/environment.yml +++ b/environment.yml @@ -61,4 +61,4 @@ dependencies: - mamba - lazy_loader - defusedxml - - neo + - python-neo