From 26a1a0180c9a6a2d1b769e0d96b2fdbe0226e9e8 Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 13:44:24 +0200 Subject: [PATCH 01/11] ENH: allow Epochs to hold trials of different duration Some experiments produce trials whose length is part of what is being measured: a gait cycle, a spoken word, a sleep stage. Cutting them to a common window either pads the short ones or truncates the rest, and both choices are made silently. `tmin` and `tmax` now accept an array with one entry per event, and `EpochsArray` accepts a list of (n_channels, n_times) arrays, deriving each epoch's `tmax` from its own length. Bounds that carry no actual variation collapse back to a single value, so nothing about the existing scalar path changes. The object reports itself through `variable_duration` and describes its trials with `durations` and `get_times(epoch)`. `times` refuses rather than inventing a shared axis, since returning the longest epoch's axis would leave `len(epochs.times) == data.shape[-1]` false while looking ordinary; `as_fixed()` returns the padded copy along with the number of epochs contributing at each sample, so the cost of padding is visible rather than implied. Reading from `Raw` gives each epoch its own length while keeping the drop bookkeeping intact. Discussed in gh-14206. --- doc/changes/dev/14210.newfeature.rst | 1 + doc/changes/names.inc | 1 + doc/conf.py | 4 +- mne/epochs.py | 618 +++++++++++++++++++-- mne/tests/test_epochs_variable_duration.py | 354 ++++++++++++ 5 files changed, 933 insertions(+), 45 deletions(-) create mode 100644 doc/changes/dev/14210.newfeature.rst create mode 100644 mne/tests/test_epochs_variable_duration.py diff --git a/doc/changes/dev/14210.newfeature.rst b/doc/changes/dev/14210.newfeature.rst new file mode 100644 index 00000000000..fc0cde80c82 --- /dev/null +++ b/doc/changes/dev/14210.newfeature.rst @@ -0,0 +1 @@ +Allow :class:`mne.Epochs` to hold trials of different duration by passing ``tmin`` and/or ``tmax`` as arrays with one entry per event, with :meth:`mne.Epochs.as_fixed` to obtain a fixed-duration copy spanning their union together with the number of epochs contributing at each time point, by `Sina Esmaeili`_. diff --git a/doc/changes/names.inc b/doc/changes/names.inc index 6856ca7e2da..6f48b057789 100644 --- a/doc/changes/names.inc +++ b/doc/changes/names.inc @@ -425,6 +425,7 @@ .. _Simon Kornblith: https://simonster.com .. _Simon M. Hofmann: https://github.com/SHEscher .. _Simon-Shlomo Poil: https://github.com/simon-shlomo +.. _Sina Esmaeili: https://github.com/snesmaeili .. _Sondre Foslien: https://github.com/sondrfos .. _Sophie Herbst: https://github.com/SophieHerbst .. _Sourav Singh: https://github.com/souravsingh diff --git a/doc/conf.py b/doc/conf.py index 487e86ab884..c492a91498b 100644 --- a/doc/conf.py +++ b/doc/conf.py @@ -779,8 +779,8 @@ def fix_sklearn_inherited_docstrings(app, what, name, obj, options, lines): ( # BaseRaw attributes are documented in Raw "py:obj", "(filename|metadata|proj|times|tmax|tmin|annotations|ch_names" - "|compensation_grade|duration|filenames|first_samp|first_time" - "|last_samp|n_times|proj|times|tmax|tmin)", + "|compensation_grade|duration|durations|filenames|first_samp|first_time" + "|last_samp|n_times|proj|times|tmax|tmin|variable_duration)", ), ] suppress_warnings = [ diff --git a/mne/epochs.py b/mne/epochs.py index f313a94e2ac..98897bf2c9f 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -385,6 +385,176 @@ def _handle_event_repeated(events, event_id, event_repeated, selection, drop_log return new_events, event_id, selection, drop_log +def _check_variable_bounds(tmin, tmax, n_events): + """Normalize ``tmin``/``tmax``, which may be given per event. + + Parameters + ---------- + tmin : float | array of float + Start time(s) in seconds. + tmax : float | array of float + End time(s) in seconds. + n_events : int + Number of events, used to check the length of array inputs. + + Returns + ------- + tmin : float | array of float + The validated start time(s). + tmax : float | array of float + The validated end time(s). + variable_duration : bool + Whether either bound was given per event. + """ + arrays = [np.ndim(bound) > 0 for bound in (tmin, tmax)] + if not any(arrays): + return tmin, tmax, False + + out = [] + for name, bound in (("tmin", tmin), ("tmax", tmax)): + bound = np.atleast_1d(np.asarray(bound, dtype=float)) + if bound.ndim != 1: + raise ValueError(f"{name} must be a scalar or 1D array, got {bound.ndim}D.") + if bound.size == 1: + bound = np.repeat(bound, n_events) + elif bound.size != n_events: + raise ValueError( + f"{name} has {bound.size} entries but there are {n_events} " + "events; per-event bounds must match the number of events." + ) + if not np.all(np.isfinite(bound)): + raise ValueError(f"{name} must be finite.") + out.append(bound) + + tmin, tmax = out + if n_events > 0 and np.all(tmin == tmin[0]) and np.all(tmax == tmax[0]): + # arrays were passed, but every epoch has the same bounds, so this is + # the ordinary fixed-duration case and should behave identically + return float(tmin[0]), float(tmax[0]), False + if np.any(tmin > tmax): + bad = int(np.argmax(tmin > tmax)) + raise ValueError( + f"tmin has to be less than or equal to tmax, but epoch {bad} has " + f"tmin={tmin[bad]:g} > tmax={tmax[bad]:g}." + ) + return tmin, tmax, True + + +def _check_variable_data(data, tmin, tmax, events, sfreq): + """Validate per-epoch data for variable-duration epochs. + + Parameters + ---------- + data : list of array + Per-epoch data, each of shape ``(n_channels, n_times_i)``. + tmin : array of float + Per-epoch start times in seconds. + tmax : array of float + Per-epoch end times in seconds. + events : array of int + The events. + sfreq : float + The sampling frequency in Hz. + """ + if isinstance(data, np.ndarray) and data.ndim == 3: + data = list(data) + if not isinstance(data, list | tuple): + raise TypeError( + "With per-event tmin/tmax, data must be a list of arrays of shape " + f"(n_channels, n_times_i), got {type(data).__name__}." + ) + if len(data) != len(events): + raise ValueError( + f"data has {len(data)} epochs but there are {len(events)} events." + ) + n_channels = None + for ii, epoch in enumerate(data): + epoch = np.asarray(epoch) + if epoch.ndim != 2: + raise ValueError( + f"Epoch {ii} has ndim={epoch.ndim}, expected 2 (n_channels, n_times)." + ) + if n_channels is None: + n_channels = epoch.shape[0] + elif epoch.shape[0] != n_channels: + raise ValueError( + f"Epoch {ii} has {epoch.shape[0]} channels but epoch 0 has " + f"{n_channels}. Only the time axis may vary between epochs." + ) + want = round((tmax[ii] - tmin[ii]) * sfreq) + 1 + if epoch.shape[1] != want: + raise ValueError( + f"Epoch {ii} has {epoch.shape[1]} samples but tmin={tmin[ii]:g} " + f"and tmax={tmax[ii]:g} at {sfreq:g} Hz imply {want}." + ) + + +def _check_variable_unsupported(*, baseline, reject_tmin, reject_tmax, decim, preload): + """Reject options this implementation does not yet handle. + + Parameters + ---------- + baseline : tuple | None + The requested baseline. + reject_tmin : float | None + Start of the rejection window. + reject_tmax : float | None + End of the rejection window. + decim : int + Decimation factor. + preload : bool + Whether the data is preloaded. + """ + if baseline is not None: + raise NotImplementedError( + "Baseline correction is not implemented for variable-duration " + "epochs; pass baseline=None and correct afterwards. The baseline " + "window would have to be checked against each epoch separately." + ) + for name, value in (("reject_tmin", reject_tmin), ("reject_tmax", reject_tmax)): + if value is not None: + raise NotImplementedError( + f"{name} is not implemented for variable-duration epochs, " + "because the window is not guaranteed to exist in every epoch." + ) + if decim != 1: + raise NotImplementedError( + "decim is not implemented for variable-duration epochs." + ) + if not preload: + raise NotImplementedError( + "Variable-duration epochs must be preloaded. Reading them lazily " + "from disk requires per-epoch bounds in the I/O layer, which is " + "not part of this change." + ) + + +def _is_variable_duration_data(data): + """Whether ``data`` is a sequence of arrays with differing lengths. + + Parameters + ---------- + data : array | list of array + Candidate epoch data. + + Returns + ------- + variable : bool + True if the entries cannot share one time axis. + """ + if isinstance(data, np.ndarray): + return False + if not isinstance(data, list | tuple) or len(data) == 0: + return False + lengths = set() + for epoch in data: + epoch = np.asarray(epoch) + if epoch.ndim != 2: + return False + lengths.add(epoch.shape[1]) + return len(lengths) > 1 + + @fill_doc class BaseEpochs( ProjMixin, @@ -410,9 +580,11 @@ class BaseEpochs( Parameters ---------- %(info_not_none)s - data : ndarray | None + data : ndarray | list of ndarray | None If ``None``, data will be read from the Raw object. If ndarray, must be - of shape (n_epochs, n_channels, n_times). + of shape (n_epochs, n_channels, n_times). A list of + (n_channels, n_times) arrays, one per epoch, gives epochs of differing + duration. %(events_epochs)s %(event_id)s %(epochs_tmin_tmax)s @@ -459,15 +631,22 @@ class BaseEpochs( used as a constructor for Epochs objects (use instead :class:`mne.Epochs`). """ + # Variable-duration epochs keep one array per epoch in this same slot. The + # declaration describes the rectangular case that every reshaping code path + # below relies on; those paths refuse before they run when durations vary. + _data: np.ndarray | None + _tmin_per_epoch: np.ndarray + _tmax_per_epoch: np.ndarray + @verbose def __init__( self, info: Info, - data: np.ndarray | None, + data: np.ndarray | list[np.ndarray] | None, events: np.ndarray, event_id: int | list[int] | dict | str | list[str] | None = None, - tmin: float = -0.2, - tmax: float = 0.5, + tmin: float | np.ndarray = -0.2, + tmax: float | np.ndarray = 0.5, baseline: tuple[float | None, float | None] | None = (None, 0), raw: "BaseRaw | list | None" = None, picks: str | np.ndarray | slice | None = None, @@ -597,6 +776,14 @@ def __init__( self.metadata = metadata # do not set self.events here, let subclass do it + # Variable-duration epochs: tmin and/or tmax may be given per event. The + # scalar path below is untouched; only _variable_duration is new. + tmin, tmax, self._variable_duration = _check_variable_bounds( + tmin, + tmax, + len(self.events) if getattr(self, "events", None) is not None else 0, + ) + if (detrend not in [None, 0, 1]) or isinstance(detrend, bool): raise ValueError("detrend must be None, 0, or 1") self.detrend = detrend @@ -614,8 +801,16 @@ def __init__( self.preload = False self._data = None self._do_baseline = True + elif self._variable_duration: + _check_variable_data(data, tmin, tmax, self.events, self.info["sfreq"]) + self.preload = True + self._data = [ # ty: ignore[invalid-assignment] # ragged payload + np.asarray(epoch, dtype=np.float64) for epoch in data + ] + self._do_baseline = False else: assert decim == 1 + assert not isinstance(data, list) # only the branch above takes a list if ( data.ndim != 3 or data.shape[2] != round((tmax - tmin) * self.info["sfreq"]) + 1 @@ -630,13 +825,35 @@ def __init__( self._do_baseline = False self._offset = None - if tmin > tmax: + if not self._variable_duration and tmin > tmax: raise ValueError("tmin has to be less than or equal to tmax") # Handle times sfreq = float(self.info["sfreq"]) - start_idx = int(round(tmin * sfreq)) - self._raw_times = np.arange(start_idx, int(round(tmax * sfreq)) + 1) / sfreq + if self._variable_duration: + # Per-epoch bounds are kept, and `times` spans their union so that + # it agrees with `as_fixed()`. Per-epoch axes are available from + # `get_times()`; see `durations`. + self._tmin_per_epoch = tmin + self._tmax_per_epoch = tmax + start_idx = int(round(tmin.min() * sfreq)) + stop_idx = int(round(tmax.max() * sfreq)) + _check_variable_unsupported( + baseline=baseline, + reject_tmin=reject_tmin, + reject_tmax=reject_tmax, + decim=decim, + preload=self.preload or preload_at_end, + ) + baseline = None + else: + # only the variable-duration path carries per-epoch bounds; every + # reader of these guards on _variable_duration first + self._tmin_per_epoch = None # ty: ignore[invalid-assignment] + self._tmax_per_epoch = None # ty: ignore[invalid-assignment] + start_idx = int(round(tmin * sfreq)) + stop_idx = int(round(tmax * sfreq)) + self._raw_times = np.arange(start_idx, stop_idx + 1) / sfreq self._set_times(self._raw_times) # check reject_tmin and reject_tmax @@ -670,12 +887,21 @@ def __init__( # decimation self._decim = 1 - self.decimate(decim) + if not self._variable_duration: + self.decimate(decim) + else: + # decim != 1 is rejected above; decimate() would densify _data + self._decim_slice = slice(None, None, None) # baseline correction: replace `None` tuple elements with actual times - self.baseline = _check_baseline( - baseline, times=self.times, sfreq=self.info["sfreq"] - ) + if self._variable_duration: + # baseline is forced to None above, and there is no shared time axis + # to resolve None endpoints against + self.baseline = None + else: + self.baseline = _check_baseline( + baseline, times=self.times, sfreq=self.info["sfreq"] + ) if self.baseline is not None and self.baseline != baseline: logger.info( f"Setting baseline interval to " @@ -743,7 +969,9 @@ def _check_consistency(self): assert len(self.drop_log) >= len(self.events) assert len(self.selection) == sum(len(dl) == 0 for dl in self.drop_log) assert hasattr(self, "_times_readonly") - assert not self.times.flags["WRITEABLE"] + # `times` is deliberately unavailable when durations vary, so check the + # underlying vector rather than going through the property + assert not self._times_readonly.flags["WRITEABLE"] assert isinstance(self.drop_log, tuple) assert all(isinstance(log, tuple) for log in self.drop_log) assert all(isinstance(s, str) for log in self.drop_log for s in log) @@ -773,6 +1001,168 @@ def reset_index(self) -> None: self.drop_log = (tuple(),) * len(self.events) self._check_consistency() + # -- variable-duration epochs --------------------------------------- + @property + def variable_duration(self): + """Whether the epochs have per-event ``tmin``/``tmax``.""" + return self._variable_duration + + @property + def times(self): + """The time axis shared by all epochs, in seconds. + + Raises + ------ + RuntimeError + When durations vary, because there is no such axis. Returning the + union would leave ``len(epochs.times) == data.shape[-1]`` false while + looking ordinary, which is how silent errors happen downstream. Use + :meth:`get_times` for one epoch, or :meth:`as_fixed` for the padded + common axis. + """ + if self._variable_duration: + raise RuntimeError( + f"These {len(self.events)} epochs have durations from " + f"{self.durations.min():.3f} to {self.durations.max():.3f} s, so " + "there is no time axis they share. Use get_times(epoch) for one " + "epoch's axis, durations for their lengths, or as_fixed() for " + "the padded common axis." + ) + return super().times + + @property + def tmin(self): + """First time point, or one per epoch when durations vary.""" + if self._variable_duration: + return self._tmin_per_epoch + return self.times[0] + + @property + def tmax(self): + """Last time point, or one per epoch when durations vary.""" + if self._variable_duration: + return self._tmax_per_epoch + return self.times[-1] + + @property + def durations(self): + """Duration of each epoch in seconds.""" + if not self._variable_duration: + return np.full(len(self.events), self.times[-1] - self.times[0]) + return self._tmax_per_epoch - self._tmin_per_epoch + + def get_times(self, epoch=None): + """Return the time vector of one epoch, or of all of them. + + Parameters + ---------- + epoch : int | None + Index of the epoch. If ``None``, return one time vector per epoch. + + Returns + ------- + times : array | list of array + Time vector(s) in seconds. For fixed-duration epochs every entry is + :attr:`times`. + """ + if epoch is None: + return [self.get_times(ii) for ii in range(len(self.events))] + if not self._variable_duration: + return self.times + sfreq = float(self.info["sfreq"]) + start = int(round(self._tmin_per_epoch[epoch] * sfreq)) + stop = int(round(self._tmax_per_epoch[epoch] * sfreq)) + return np.arange(start, stop + 1) / sfreq + + def as_fixed(self, pad_value=np.nan): + """Return fixed-duration epochs spanning the union of all epochs. + + Parameters + ---------- + pad_value : float + Value used where an epoch does not extend across the full window. + + Returns + ------- + epochs : instance of EpochsArray + Fixed-duration epochs from ``min(tmin)`` to ``max(tmax)``. + n_contributing : array, shape (n_times,) + Number of epochs carrying real data at each time point. + + Notes + ----- + Shorter epochs are padded to reach the common window, so the number of + epochs contributing to any reduction varies across the window. That + makes the noise level time-dependent, and the scalar ``nave`` carried by + :class:`mne.Evoked` cannot express it, which matters wherever ``nave`` + scales a noise covariance. ``n_contributing`` is returned rather than + discarded so this stays visible; see + :gh:`14206` for the discussion. + """ + if not self._variable_duration: + return self.copy(), np.full(len(self.times), len(self.events)) + # the union vector, which `times` deliberately refuses to hand out + times = self._times_readonly + n_times = len(times) + data = np.full( + (len(self.events), len(self.ch_names), n_times), pad_value, dtype=float + ) + n_contributing = np.zeros(n_times, int) + sfreq = float(self.info["sfreq"]) + offset = int(round(times[0] * sfreq)) + ragged = self._data + assert ragged is not None # variable-duration epochs are always preloaded + for ii, epoch in enumerate(ragged): + start = int(round(self._tmin_per_epoch[ii] * sfreq)) - offset + stop = start + epoch.shape[1] + data[ii, :, start:stop] = epoch + n_contributing[start:stop] += 1 + out = EpochsArray( + data, + self.info.copy(), + events=self.events.copy(), + event_id=self.event_id, + tmin=float(times[0]), + metadata=self.metadata, + selection=self.selection, + drop_log=self.drop_log, + verbose=False, + ) + return out, n_contributing + + def _get_variable_data(self, *, picks=None, item=None, copy=True): + """Return per-epoch data when durations vary. + + Parameters + ---------- + picks : str | array-like | slice | None + Channels to include. + item : slice | array-like | str | list | None + Epochs to include. + copy : bool + Whether to copy the data. + + Returns + ------- + data : list of array + One ``(n_channels, n_times_i)`` array per epoch. A single array is + not returned, because there is no length every epoch shares; use + :meth:`as_fixed` to obtain one. + """ + if item is None: + item = slice(None) + sel = np.arange(len(self.events))[item] if not isinstance(item, str) else None + if sel is None: + raise NotImplementedError( + "Selecting variable-duration epochs by condition name is not " + "implemented; index with integers or a slice." + ) + ch_idx = _picks_to_idx(self.info, picks, "all", exclude=()) + ragged = self._data + assert ragged is not None # variable-duration epochs are always preloaded + out = [ragged[ii][ch_idx] for ii in sel] + return [epoch.copy() for epoch in out] if copy else out + def load_data(self) -> Self: """Load the data if not already preloaded. @@ -789,13 +1179,20 @@ def load_data(self) -> Self: """ if self.preload: return self - self._data = self._get_data() + if self._variable_duration: + # a list, one array per epoch; see the note on _data above + self._data = self._load_variable_from_raw() + else: + self._data = self._get_data() self.preload = True self._do_baseline = False self._decim_slice = slice(None, None, None) self._decim = 1 - self._raw_times = self.times - assert self._data.shape[-1] == len(self.times) + if not self._variable_duration: + # variable-duration epochs already set _raw_times to the union span + # in __init__, and self.times refuses to answer for them + self._raw_times = self.times + assert self._data.shape[-1] == len(self.times) self._raw = None # shouldn't need it anymore return self @@ -953,14 +1350,81 @@ def _reject_setup(self, reject, flat, *, allow_callable=False): reject_imax = idxs[-1] self._reject_time = slice(reject_imin, reject_imax) + def _load_variable_from_raw(self): + """Read every epoch at its own length, dropping bad ones. + + Returns + ------- + data : list of array + One ``(n_channels, n_times_i)`` array per retained epoch. + + Notes + ----- + Mirrors the drop bookkeeping of :meth:`drop_bad` for the fixed-duration + path, but collects a list because there is no length the epochs share. + """ + detrend_picks = self._detrend_picks + drop_log = list(self.drop_log) + good_idx, out = [], [] + for idx, sel in enumerate(self.selection): + epoch_noproj = self._get_epoch_from_raw(idx) + epoch_noproj = self._detrend_offset_decim(epoch_noproj, detrend_picks) + epoch = self._project_epoch(epoch_noproj) + epoch_out = epoch_noproj if self._do_delayed_proj else epoch + is_good, bad_tuple = self._is_good_epoch( + epoch, n_times=self._n_times_per_epoch(idx) + ) + if not is_good: + drop_log[sel] = drop_log[sel] + bad_tuple + continue + good_idx.append(idx) + out.append(epoch_out) + + good_idx = np.asarray(good_idx, dtype=int) + if len(good_idx) != len(self.events): + self._tmin_per_epoch = self._tmin_per_epoch[good_idx] + self._tmax_per_epoch = self._tmax_per_epoch[good_idx] + self.events = self.events[good_idx] + self.selection = self.selection[good_idx] + if self.metadata is not None: + GetEpochsMixin.metadata.fset( + self, self.metadata.iloc[good_idx], verbose=False + ) + self.drop_log = tuple(drop_log) + self._bad_dropped = True + logger.info(f"{len(out)} matching events found after rejection") + return out + + def _n_times_per_epoch(self, idx): + """Return the number of samples in one epoch. + + Parameters + ---------- + idx : int + Index of the epoch. + + Returns + ------- + n_times : int + Number of samples. + """ + if not self._variable_duration: + return len(self.times) + sfreq = float(self.info["sfreq"]) + return ( + int(round((self._tmax_per_epoch[idx] - self._tmin_per_epoch[idx]) * sfreq)) + + 1 + ) + @verbose # verbose is used by mne-realtime - def _is_good_epoch(self, data, verbose=None): + def _is_good_epoch(self, data, verbose=None, *, n_times=None): """Determine if epoch is good.""" if isinstance(data, str): return False, (data,) if data is None: return False, ("NO_DATA",) - n_times = len(self.times) + if n_times is None: + n_times = len(self.times) if data.shape[1] < n_times: # epoch is too short ie at the end of the data return False, ("TOO_SHORT",) @@ -1873,7 +2337,15 @@ def _get_data( epoch = self._project_epoch(epoch_noproj) epoch_out = epoch_noproj if self._do_delayed_proj else epoch - is_good, bad_tuple = self._is_good_epoch(epoch, verbose=verbose) + is_good, bad_tuple = self._is_good_epoch( + epoch, + verbose=verbose, + n_times=( + self._n_times_per_epoch(idx) + if self._variable_duration + else None + ), + ) if not is_good: assert isinstance(bad_tuple, tuple) assert all(isinstance(x, str) for x in bad_tuple) @@ -2043,6 +2515,8 @@ def get_data( The epochs data. Will be a copy when ``copy=True`` and will be a view when possible when ``copy=False``. """ + if self._variable_duration: + return self._get_variable_data(picks=picks, item=item, copy=copy) return self._get_data( picks=picks, item=item, units=units, tmin=tmin, tmax=tmax, copy=copy ) @@ -2158,8 +2632,13 @@ def __repr__(self): """Build string representation.""" s = f"{len(self.events)} events " s += "(all good)" if self._bad_dropped else "(good & bad)" - s += f", {self.tmin:.3f}".rstrip("0").rstrip(".") - s += f" – {self.tmax:.3f}".rstrip("0").rstrip(".") + if self._variable_duration: + s += f", {self.tmin.min():.3f}".rstrip("0").rstrip(".") + s += f" – {self.tmax.max():.3f}".rstrip("0").rstrip(".") + s += f" s, durations {self.durations.min():.3f}–{self.durations.max():.3f}" + else: + s += f", {self.tmin:.3f}".rstrip("0").rstrip(".") + s += f" – {self.tmax:.3f}".rstrip("0").rstrip(".") s += " s (baseline " if self.baseline is None: s += "off" @@ -3685,8 +4164,8 @@ def __init__( raw: "BaseRaw", events: np.ndarray | None = None, event_id: int | list[int] | dict | str | list[str] | None = None, - tmin: float = -0.2, - tmax: float = 0.5, + tmin: float | np.ndarray = -0.2, + tmax: float | np.ndarray = 0.5, baseline: tuple[float | None, float | None] | None = (None, 0), picks: str | np.ndarray | slice | None = None, preload: bool = False, @@ -3789,24 +4268,37 @@ def _get_epoch_from_raw(self, idx, verbose=None): assert not isinstance(self._raw, list) # a single Raw, unlike EpochsFIF sfreq = self._raw.info["sfreq"] event_samp = self.events[idx, 0] - # Read a data segment from "start" to "stop" in samples + # Read a data segment from "start" to "stop" in samples. With per-event + # bounds each epoch has its own window; otherwise all epochs share + # self._raw_times. first_samp = self._raw.first_samp - start = int(round(event_samp + self._raw_times[0] * sfreq)) + if self._variable_duration: + epoch_tmin = self._tmin_per_epoch[idx] + n_samples = self._n_times_per_epoch(idx) + else: + epoch_tmin = self._raw_times[0] + n_samples = len(self._raw_times) + start = int(round(event_samp + epoch_tmin * sfreq)) start -= first_samp - stop = start + len(self._raw_times) + stop = start + n_samples # reject_tmin, and reject_tmax need to be converted to samples to # check the reject_by_annotation boundaries: reject_start, reject_stop reject_tmin = self.reject_tmin if reject_tmin is None: - reject_tmin = self._raw_times[0] + reject_tmin = epoch_tmin reject_start = int(round(event_samp + reject_tmin * sfreq)) reject_start -= first_samp + epoch_tmax = ( + self._tmax_per_epoch[idx] + if self._variable_duration + else self._raw_times[-1] + ) reject_tmax = self.reject_tmax if reject_tmax is None: - reject_tmax = self._raw_times[-1] - diff = int(round((self._raw_times[-1] - reject_tmax) * sfreq)) + reject_tmax = epoch_tmax + diff = int(round((epoch_tmax - reject_tmax) * sfreq)) reject_stop = stop - diff logger.debug(f" Getting epoch for {start}-{stop}") @@ -3881,10 +4373,10 @@ class EpochsArray(BaseEpochs): @verbose def __init__( self, - data: np.ndarray, + data: np.ndarray | list[np.ndarray], info: Info, events: np.ndarray | None = None, - tmin: float = 0.0, + tmin: float | np.ndarray = 0.0, event_id: int | list[int] | dict | str | list[str] | None = None, reject: dict | None = None, flat: dict | None = None, @@ -3901,20 +4393,60 @@ def __init__( raw_sfreq: float | None = None, verbose: bool | str | int | None = None, ): - dtype = np.complex128 if np.any(np.iscomplex(data)) else np.float64 - data = np.asanyarray(data, dtype=dtype) - if data.ndim != 3: - raise ValueError( - "Data must be a 3D array of shape (n_epochs, n_channels, n_samples)" + # A list of differently-shaped arrays means variable-duration epochs; + # each epoch's tmax follows from its own length. A ragged list cannot + # become a 3D array, so this must be decided before the cast below. + # Match the decision BaseEpochs will make: bounds that are given as + # arrays but carry no actual variation are the ordinary fixed case. + _tmin_arr = np.asarray(tmin, dtype=float) + variable = _is_variable_duration_data(data) or ( + _tmin_arr.ndim > 0 + and _tmin_arr.size > 0 + and not np.all(_tmin_arr == _tmin_arr.flat[0]) + ) + if not variable and _tmin_arr.ndim > 0: + tmin = float(_tmin_arr.flat[0]) if _tmin_arr.size else 0.0 + if variable: + data = [np.asanyarray(epoch) for epoch in data] + dtype = ( + np.complex128 + if any(np.any(np.iscomplex(epoch)) for epoch in data) + else np.float64 + ) + data = [epoch.astype(dtype, copy=False) for epoch in data] + for ii, epoch in enumerate(data): + if epoch.ndim != 2: + raise ValueError( + f"Epoch {ii} must be 2D (n_channels, n_samples), got " + f"{epoch.ndim}D." + ) + if len(info["ch_names"]) != epoch.shape[0]: + raise ValueError("Info and data must have same number of channels.") + if events is None: + events = _gen_events(len(data)) + info = info.copy() # do not modify original info + tmin = np.broadcast_to(np.asarray(tmin, dtype=float), (len(data),)).copy() + tmax = np.array( + [ + (epoch.shape[1] - 1) / info["sfreq"] + start + for epoch, start in zip(data, tmin) + ] ) + else: + dtype = np.complex128 if np.any(np.iscomplex(data)) else np.float64 + data = np.asanyarray(data, dtype=dtype) + if data.ndim != 3: + raise ValueError( + "Data must be a 3D array of shape (n_epochs, n_channels, n_samples)" + ) - if len(info["ch_names"]) != data.shape[1]: - raise ValueError("Info and data must have same number of channels.") - if events is None: - n_epochs = len(data) - events = _gen_events(n_epochs) - info = info.copy() # do not modify original info - tmax = (data.shape[2] - 1) / info["sfreq"] + tmin + if len(info["ch_names"]) != data.shape[1]: + raise ValueError("Info and data must have same number of channels.") + if events is None: + n_epochs = len(data) + events = _gen_events(n_epochs) + info = info.copy() # do not modify original info + tmax = (data.shape[2] - 1) / info["sfreq"] + tmin super().__init__( info, diff --git a/mne/tests/test_epochs_variable_duration.py b/mne/tests/test_epochs_variable_duration.py new file mode 100644 index 00000000000..04d0e559e0c --- /dev/null +++ b/mne/tests/test_epochs_variable_duration.py @@ -0,0 +1,354 @@ +"""Tests for epochs whose trials have different durations.""" + +# Authors: The MNE-Python contributors. +# License: BSD-3-Clause +# Copyright the MNE-Python contributors. + +import numpy as np +import pytest +from numpy.testing import assert_allclose, assert_array_equal + +from mne import EpochsArray, create_info + +SFREQ = 100.0 +CH_NAMES = ["a", "b", "c"] + + +def _events(n_epochs): + """Return a simple event array.""" + return np.c_[ + np.arange(n_epochs) * 500 + 200, + np.zeros(n_epochs, int), + np.ones(n_epochs, int), + ] + + +def _make(tmin, tmax, seed=0): + """Build variable-duration epochs from per-event bounds.""" + rng = np.random.default_rng(seed) + tmin = np.asarray(tmin, dtype=float) + tmax = np.asarray(tmax, dtype=float) + data = [ + rng.standard_normal((len(CH_NAMES), round((b - a) * SFREQ) + 1)) * 1e-6 + for a, b in zip(tmin, tmax) + ] + info = create_info(CH_NAMES, SFREQ, "eeg") + return EpochsArray( + data, + info, + events=_events(len(data)), + tmin=tmin, + baseline=None, + verbose=False, + ) + + +@pytest.fixture +def variable(): + """Four epochs sharing tmin with different tmax.""" + return _make(np.full(4, -0.2), [0.5, 0.9, 0.7, 0.6]) + + +# -- the scalar path must be untouched ------------------------------------ +def test_scalar_path_unchanged(): + """Test that fixed-duration epochs behave exactly as before.""" + rng = np.random.default_rng(0) + data = rng.standard_normal((5, len(CH_NAMES), 71)) * 1e-6 + info = create_info(CH_NAMES, SFREQ, "eeg") + epochs = EpochsArray(data, info, events=_events(5), tmin=-0.2, verbose=False) + + assert not epochs.variable_duration + assert isinstance(epochs.tmin, float) + assert isinstance(epochs.tmax, float) + assert_allclose(epochs.tmin, -0.2) + assert epochs.get_data().shape == (5, len(CH_NAMES), 71) + assert epochs.average().data.shape == (len(CH_NAMES), 71) + # durations is new but must be defined for the fixed case too + assert_allclose(epochs.durations, np.full(5, 0.7)) + + +def test_equal_bounds_arrays_match_scalar(): + """Test that per-event bounds that happen to be equal match the scalar path.""" + rng = np.random.default_rng(3) + data = rng.standard_normal((4, len(CH_NAMES), 71)) * 1e-6 + info = create_info(CH_NAMES, SFREQ, "eeg") + + scalar = EpochsArray(data, info, events=_events(4), tmin=-0.2, verbose=False) + arrays = EpochsArray( + list(data), info, events=_events(4), tmin=np.full(4, -0.2), verbose=False + ) + + # a list of equal-length arrays is not ragged, so this stays the scalar path + assert not arrays.variable_duration + assert_allclose(arrays.times, scalar.times) + assert_allclose(arrays.get_data(), scalar.get_data()) + + +# -- construction --------------------------------------------------------- +def test_bounds_and_durations(variable): + """Test the per-epoch bounds derived from the data.""" + assert variable.variable_duration + assert_allclose(variable.tmin, np.full(4, -0.2)) + assert_allclose(variable.tmax, [0.5, 0.9, 0.7, 0.6], atol=1e-12) + assert_allclose(variable.durations, [0.7, 1.1, 0.9, 0.8], atol=1e-12) + + +def test_per_epoch_time_vectors(variable): + """Test that each epoch reports its own axis, since none is shared.""" + per_epoch = variable.get_times() + assert len(per_epoch) == 4 + for ii, epoch_times in enumerate(per_epoch): + assert_allclose(epoch_times[0], variable.tmin[ii]) + assert_allclose(epoch_times[-1], variable.tmax[ii], atol=1e-12) + assert len(epoch_times) == variable.get_data()[ii].shape[-1] + + +def test_get_data_returns_one_array_per_epoch(variable): + """Test that data comes back per epoch, since no common length exists.""" + data = variable.get_data() + assert isinstance(data, list) + lengths = [epoch.shape[1] for epoch in data] + assert lengths == [71, 111, 91, 81] + assert all(epoch.shape[0] == len(CH_NAMES) for epoch in data) + + picked = variable.get_data(picks=["a", "c"]) + assert all(epoch.shape[0] == 2 for epoch in picked) + + +def test_only_time_may_vary(): + """Test that a varying channel count is rejected.""" + info = create_info(CH_NAMES, SFREQ, "eeg") + data = [np.zeros((3, 71)), np.zeros((2, 111))] + with pytest.raises(ValueError, match="same number of channels"): + EpochsArray(data, info, events=_events(2), tmin=np.zeros(2), verbose=False) + + +def test_bounds_validation(): + """Test the checks on per-event bounds.""" + from mne.epochs import _check_variable_bounds + + tmin, tmax, variable = _check_variable_bounds(-0.2, 0.5, 3) + assert not variable + assert tmin == -0.2 + + tmin, tmax, variable = _check_variable_bounds(-0.2, np.array([0.5, 0.9, 0.7]), 3) + assert variable + assert_array_equal(tmin, np.full(3, -0.2)) + + with pytest.raises(ValueError, match="entries but there are"): + _check_variable_bounds(np.zeros(2), 0.5, 3) + with pytest.raises(ValueError, match="less than or equal to tmax"): + _check_variable_bounds(np.array([0.5, 0.0]), np.array([0.1, 1.0]), 2) + with pytest.raises(ValueError, match="must be finite"): + _check_variable_bounds(np.array([np.nan, 0.0]), np.array([1.0, 1.0]), 2) + + +@pytest.mark.parametrize( + "kwargs, match", + [ + (dict(baseline=(None, 0)), "Baseline correction is not implemented"), + (dict(reject_tmin=-0.1), "reject_tmin is not implemented"), + (dict(reject_tmax=0.1), "reject_tmax is not implemented"), + ], +) +def test_unsupported_options_raise(kwargs, match): + """Test that options outside this implementation are refused clearly.""" + rng = np.random.default_rng(0) + tmin = np.full(2, -0.2) + data = [rng.standard_normal((len(CH_NAMES), n)) for n in (71, 111)] + info = create_info(CH_NAMES, SFREQ, "eeg") + with pytest.raises(NotImplementedError, match=match): + EpochsArray(data, info, events=_events(2), tmin=tmin, verbose=False, **kwargs) + + +# -- as_fixed ------------------------------------------------------------- +def test_as_fixed_spans_union_and_reports_support(variable): + """Test padding to a common window and the contributor count it implies.""" + fixed, n_contributing = variable.as_fixed() + + assert not fixed.variable_duration + n_union = max(epoch.shape[-1] for epoch in variable.get_data()) + assert fixed.get_data().shape == (4, len(CH_NAMES), n_union) + assert_allclose(fixed.times[0], variable.tmin.min()) + assert_allclose(fixed.times[-1], variable.tmax.max()) + + # every epoch covers the start; only the longest reaches the end + assert n_contributing[0] == 4 + assert n_contributing[-1] == 1 + assert n_contributing.min() < n_contributing.max() + # support never grows once epochs start dropping out + peak = int(np.argmax(n_contributing)) + assert np.all(np.diff(n_contributing[peak:]) <= 0) + + +def test_as_fixed_preserves_data_and_pads_the_rest(variable): + """Test that real samples survive and the remainder is marked.""" + fixed, _ = variable.as_fixed() + dense = fixed.get_data() + for ii, epoch in enumerate(variable.get_data()): + n = epoch.shape[1] + assert_allclose(dense[ii, :, :n], epoch) + assert np.isnan(dense[ii, :, n:]).all() + + +def test_as_fixed_on_fixed_epochs_is_a_copy(): + """Test that as_fixed is defined, and trivial, for fixed-duration epochs.""" + rng = np.random.default_rng(1) + data = rng.standard_normal((3, len(CH_NAMES), 71)) * 1e-6 + info = create_info(CH_NAMES, SFREQ, "eeg") + epochs = EpochsArray(data, info, events=_events(3), tmin=-0.2, verbose=False) + + fixed, n_contributing = epochs.as_fixed() + assert_allclose(fixed.get_data(), epochs.get_data()) + assert_array_equal(n_contributing, np.full(len(epochs.times), 3)) + + +# -- dispatch -------------------------------------------------------------- + + +# -- operations that stay native ------------------------------------------- + + +# -- the time axis --------------------------------------------------------- +def test_times_refuses_to_invent_a_shared_axis(variable): + """Test that ``times`` raises rather than handing back the union. + + ``len(epochs.times) == data.shape[-1]`` has always held. Returning the union + would leave it false while looking ordinary. + """ + with pytest.raises(RuntimeError, match="no time axis they share"): + variable.times + with pytest.raises(RuntimeError, match="get_times"): + variable.times + + +def test_union_axis_is_reachable_through_as_fixed(variable): + """Test that the padded common axis is available when asked for.""" + fixed, _ = variable.as_fixed() + assert_allclose(fixed.times[0], variable.tmin.min()) + assert_allclose(fixed.times[-1], variable.tmax.max()) + assert len(fixed.times) == fixed.get_data().shape[-1] + + +def test_fixed_epochs_still_have_times(): + """Test that the scalar path is untouched by any of this.""" + rng = np.random.default_rng(2) + info = create_info(CH_NAMES, SFREQ, "eeg") + epochs = EpochsArray( + rng.standard_normal((3, len(CH_NAMES), 71)) * 1e-6, + info, + events=_events(3), + tmin=-0.2, + verbose=False, + ) + assert len(epochs.times) == 71 + assert epochs.average().data.shape == (len(CH_NAMES), 71) + + +# -- construction from Raw ------------------------------------------------- +def _raw(n_seconds=30.0, sfreq=SFREQ): + """Return a small continuous recording.""" + from mne.io import RawArray + + rng = np.random.default_rng(3) + info = create_info(CH_NAMES, sfreq, "eeg") + return RawArray( + rng.standard_normal((len(CH_NAMES), int(n_seconds * sfreq))) * 1e-6, + info, + verbose=False, + ) + + +def test_from_raw_preserves_each_slice(): + """Test that epochs read from Raw match the samples they came from.""" + from mne import Epochs + + raw = _raw() + onsets = np.array([100, 700, 1400, 2100]) + events = np.c_[onsets, np.zeros(4, int), np.ones(4, int)] + tmax = np.array([0.5, 1.9, 0.9, 2.4]) + + epochs = Epochs( + raw, + events, + {"a": 1}, + tmin=np.zeros(4), + tmax=tmax, + baseline=None, + preload=True, + verbose=False, + ) + + assert epochs.variable_duration + assert len(epochs) == 4 + # +1 because both endpoints are inclusive + want_n = np.round(tmax * SFREQ).astype(int) + 1 + data = epochs.get_data() + assert [d.shape[-1] for d in data] == list(want_n) + + raw_data = raw.get_data() + for i, start in enumerate(onsets): + assert_array_equal(data[i], raw_data[:, start : start + want_n[i]]) + + +def test_from_raw_load_data_leaves_raw_times_alone(): + """Test that loading from Raw does not ask a ragged object for one axis.""" + from mne import Epochs + + raw = _raw() + events = np.c_[[100, 900], [0, 0], [1, 1]] + # load_data() used to raise here from `self._raw_times = self.times` + epochs = Epochs( + raw, + events, + {"a": 1}, + tmin=np.zeros(2), + tmax=np.array([0.5, 2.0]), + baseline=None, + preload=True, + verbose=False, + ) + assert epochs.preload + # _raw_times spans the union, which is what as_fixed() lays epochs onto + fixed, _ = epochs.as_fixed() + assert len(epochs._raw_times) == fixed.get_data().shape[-1] + + +def test_from_raw_requires_preload(): + """Test that lazy reading is refused with a reason, not an internal error.""" + from mne import Epochs + + raw = _raw() + events = np.c_[[100, 900], [0, 0], [1, 1]] + with pytest.raises(NotImplementedError, match="must be preloaded"): + Epochs( + raw, + events, + {"a": 1}, + tmin=np.zeros(2), + tmax=np.array([0.5, 2.0]), + baseline=None, + preload=False, + verbose=False, + ) + + +def test_from_raw_scalar_bounds_still_scalar(): + """Test that equal per-event bounds do not switch on the ragged path.""" + from mne import Epochs + + raw = _raw() + events = np.c_[[100, 900], [0, 0], [1, 1]] + epochs = Epochs( + raw, + events, + {"a": 1}, + tmin=np.zeros(2), + tmax=np.full(2, 0.5), + baseline=None, + preload=True, + verbose=False, + ) + assert not epochs.variable_duration + assert isinstance(epochs.tmin, float) + assert epochs.get_data().shape == (2, len(CH_NAMES), 51) From 34809b98a2ba58e575d7477e82291cf67f14b719 Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 13:45:51 +0200 Subject: [PATCH 02/11] ENH: support safe operations on variable-duration Epochs With trials of differing length the methods divide into three kinds, and guessing which one you are calling is how a wrong answer gets returned quietly. Selecting epochs, selecting channels, dropping and shifting the time origin do not care how long each trial is, so they work as they always did; the per-epoch bounds travel with the epochs they describe. That needs one branch each in GetEpochsMixin._getitem, shift_time and _pick_drop_channels, since those hold the data as one array. _pick_drop_channels replaces the list contents rather than the attribute, which keeps `_data` an ndarray for Raw, Evoked and the rest. Reductions across a shared time axis decline and say what they would need: padding makes the number of contributing epochs a function of time, which no single nave describes. Measuring it is what settled this - on 43 epochs spanning 2.0-3.6 s, average() returns an Evoked that is 44% NaN while nave reports 43 where 3 epochs remain. Per-trial operations with no ragged implementation decline too, rather than running on a padded copy and returning a wrong answer instead of a slow one. plot() is among them for now; the next commit implements it. to_data_frame keeps a warning fallback, since its result is only read. --- mne/channels/channels.py | 9 +- mne/epochs.py | 161 ++++++++++++++++++++- mne/tests/test_epochs_variable_duration.py | 149 +++++++++++++++++++ mne/utils/mixin.py | 39 ++++- 4 files changed, 352 insertions(+), 6 deletions(-) diff --git a/mne/channels/channels.py b/mne/channels/channels.py index d0e9ae73d6e..b4e40ad2b71 100644 --- a/mne/channels/channels.py +++ b/mne/channels/channels.py @@ -636,7 +636,14 @@ def _pick_drop_channels(self, idx, *, verbose=None): else: # All others (Evoked, Epochs, Raw) have chs axis=-2 axis = -2 if hasattr(self, "_data"): # skip non-preloaded Raw - self._data = self._data.take(idx, axis=axis) + if isinstance(self._data, list): + # variable-duration epochs: one array per epoch, channels are + # regular within each, so the pick applies the same way to all. + # Replacing the contents rather than the attribute keeps `_data` + # an ndarray everywhere else this mixin is used. + self._data[:] = [epoch.take(idx, axis=axis) for epoch in self._data] + else: + self._data = self._data.take(idx, axis=axis) else: assert isinstance(self, BaseRaw) and not self.preload diff --git a/mne/epochs.py b/mne/epochs.py index 98897bf2c9f..a5383939fa5 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -10,7 +10,7 @@ from collections import Counter from collections.abc import Callable, Iterable, Iterator from copy import deepcopy -from functools import partial +from functools import partial, wraps from inspect import getfullargspec from pathlib import Path from typing import TYPE_CHECKING, Literal @@ -4042,6 +4042,165 @@ def _events_from_annotations(raw, events, event_id, annotations, on_missing): return events, event_id, annotations +#: Methods whose result is only looked at. They warn and run on ``as_fixed()``, +#: which is enough for inspection because the padding is visible to whoever is +#: looking. Anything numeric is not in this table: see the two below. +_VARIABLE_FALLBACK = { + "to_data_frame": "", +} + +#: Methods that combine epochs across a shared time axis. Padding them makes the +#: number of contributing epochs a function of time, which no scalar ``nave`` can +#: describe, so they ask for a policy instead of inventing one. See :gh:`14206`. +_VARIABLE_NEEDS_POLICY = { + "average": "averaging", + "standard_error": "estimating the standard error", + "subtract_evoked": "subtracting an evoked response", + "iter_evoked": "iterating as evoked responses", + "compute_tfr": "computing a time-frequency representation", + "compute_psd": "computing a spectrum", +} + +#: Methods that are mathematically per-trial and simply have no ragged +#: implementation yet. Running them on a padded copy would return a wrong answer +#: rather than a slow one, so they raise until implemented natively. +_VARIABLE_NOT_IMPLEMENTED = { + "filter": "filtering", + "plot": "browsing", + "apply_function": "applying a function", + "apply_baseline": "baseline correction", + "crop": "cropping", + "decimate": "decimation", + "resample": "resampling", + "save": "writing to FIF", + "export": "exporting", + # these render an image over one axis and cannot draw the NaN padding that + # as_fixed() introduces, so the fallback has nothing useful to show + "plot_image": "plotting as an image", + "plot_topo_image": "plotting as a topographic image", +} + + +def _wrap_variable_fallback(func, name, note): + """Warn and fall back to ``as_fixed()`` for variable-duration epochs. + + Parameters + ---------- + func : callable + The original method. + name : str + Its name, used to look it up on the fixed-duration copy. + note : str + Extra sentence appended to the warning, or an empty string. + + Returns + ------- + wrapper : callable + The wrapped method. + """ + + @wraps(func) + def wrapper(self, *args, **kwargs): + if not getattr(self, "_variable_duration", False): + return func(self, *args, **kwargs) + message = ( + f"{name}() needs one time axis, which these variable-duration " + f"epochs do not have, so it ran on as_fixed(): every epoch padded " + f"to span {self.tmin.min():g} to {self.tmax.max():g} s." + ) + if note: + message += " " + note + message += " Call as_fixed() yourself to make this explicit." + warn(message, RuntimeWarning) + fixed, _ = self.as_fixed() + return getattr(fixed, name)(*args, **kwargs) + + return wrapper + + +def _raise_needs_policy(func, name, what): + """Raise for reductions across a time axis the epochs do not share. + + Parameters + ---------- + func : callable + The original method. + name : str + Its name. + what : str + Short description of the operation, used in the message. + + Returns + ------- + wrapper : callable + The wrapped method. + """ + + @wraps(func) + def wrapper(self, *args, **kwargs): + if not getattr(self, "_variable_duration", False): + return func(self, *args, **kwargs) + raise NotImplementedError( + f"{name}() combines epochs across a shared time axis, and these " + f"epochs do not share one. {what.capitalize()} them needs an " + "explicit policy, because the number of contributing epochs varies " + "across the window and no single nave describes it. Either call " + "as_fixed(), which pads to the union window and returns that count " + "alongside the data, or align the epochs first. See " + "https://github.com/mne-tools/mne-python/issues/14206." + ) + + return wrapper + + +def _raise_not_implemented(func, name, what): + """Raise for per-trial operations with no ragged implementation yet. + + Parameters + ---------- + func : callable + The original method. + name : str + Its name. + what : str + Short description of the operation, used in the message. + + Returns + ------- + wrapper : callable + The wrapped method. + """ + + @wraps(func) + def wrapper(self, *args, **kwargs): + if not getattr(self, "_variable_duration", False): + return func(self, *args, **kwargs) + raise NotImplementedError( + f"{name}() is not implemented for variable-duration epochs. " + f"{what.capitalize()} is per-trial and could work here, but running " + "it on a padded copy would change the result rather than just slow " + "it down, so it raises until implemented. See " + "https://github.com/mne-tools/mne-python/issues/14206." + ) + + return wrapper + + +for _name, _note in _VARIABLE_FALLBACK.items(): + _orig = getattr(BaseEpochs, _name, None) + if _orig is not None: + setattr(BaseEpochs, _name, _wrap_variable_fallback(_orig, _name, _note)) +for _name, _what in _VARIABLE_NEEDS_POLICY.items(): + _orig = getattr(BaseEpochs, _name, None) + if _orig is not None: + setattr(BaseEpochs, _name, _raise_needs_policy(_orig, _name, _what)) +for _name, _what in _VARIABLE_NOT_IMPLEMENTED.items(): + _orig = getattr(BaseEpochs, _name, None) + if _orig is not None: + setattr(BaseEpochs, _name, _raise_not_implemented(_orig, _name, _what)) +del _name, _note, _what, _orig + + @fill_doc class Epochs(BaseEpochs): """Epochs extracted from a Raw instance. diff --git a/mne/tests/test_epochs_variable_duration.py b/mne/tests/test_epochs_variable_duration.py index 04d0e559e0c..a13358f9ca9 100644 --- a/mne/tests/test_epochs_variable_duration.py +++ b/mne/tests/test_epochs_variable_duration.py @@ -9,6 +9,11 @@ from numpy.testing import assert_allclose, assert_array_equal from mne import EpochsArray, create_info +from mne.epochs import ( + _VARIABLE_FALLBACK, + _VARIABLE_NEEDS_POLICY, + _VARIABLE_NOT_IMPLEMENTED, +) SFREQ = 100.0 CH_NAMES = ["a", "b", "c"] @@ -204,9 +209,91 @@ def test_as_fixed_on_fixed_epochs_is_a_copy(): # -- dispatch -------------------------------------------------------------- +@pytest.mark.parametrize("meth", sorted(_VARIABLE_NEEDS_POLICY)) +def test_reductions_ask_for_a_policy(variable, meth): + """Test that combining epochs across a time axis they lack is refused. + + Padding first and reducing afterwards is not a slower answer, it is a + different one: one short epoch turns a whole time point into NaN, and the + scalar ``nave`` keeps reporting the full count. + """ + with pytest.raises(NotImplementedError, match="explicit policy"): + result = getattr(variable, meth)() + list(result) # iter_evoked is a generator + + +def test_policy_message_names_the_varying_count(variable): + """Test that the refusal explains itself rather than just declining.""" + with pytest.raises(NotImplementedError) as excinfo: + variable.average() + message = str(excinfo.value) + assert "varies across the window" in message + assert "as_fixed" in message + + +def test_compute_tfr_does_not_silently_pad(variable): + """Test that the transform is not quietly given padded data. + + Padding before a time-frequency transform is the opposite of the order this + work argues for, which is to transform at native duration and warp the + result. Doing it silently inside ``compute_tfr`` would ship the thing being + argued against. + """ + with pytest.raises(NotImplementedError, match="explicit policy"): + variable.compute_tfr("morlet", freqs=np.arange(10.0, 20.0, 2.0), n_cycles=2) + + +@pytest.mark.parametrize("meth", sorted(_VARIABLE_NOT_IMPLEMENTED)) +def test_per_trial_methods_raise_until_implemented(variable, meth): + """Test that per-trial work refuses rather than running on a padded copy.""" + with pytest.raises(NotImplementedError, match="not implemented"): + getattr(variable, meth)() + + +@pytest.mark.parametrize("meth", sorted(_VARIABLE_FALLBACK)) +def test_display_methods_warn_and_fall_back(variable, meth): + """Test that the remaining inspection method degrades rather than refuses.""" + if meth == "to_data_frame": + pytest.importorskip("pandas") + with pytest.warns(RuntimeWarning, match="ran on as_fixed"): + assert getattr(variable, meth)() is not None # -- operations that stay native ------------------------------------------- +def test_pick_keeps_durations(variable): + """Test that channel selection leaves the time axis alone.""" + before = variable.durations.copy() + picked = variable.copy().pick(["a", "c"]) + assert picked.ch_names == ["a", "c"] + assert_allclose(picked.durations, before) + for epoch in picked.get_data(): + assert epoch.shape[0] == 2 + + +def test_getitem_keeps_per_epoch_bounds(variable): + """Test that indexing carries the bounds with the epochs.""" + subset = variable[[0, 2]] + assert len(subset) == 2 + assert_allclose(subset.durations, variable.durations[[0, 2]]) + for got, want in zip(subset.get_data(), [variable.get_data()[i] for i in (0, 2)]): + assert_array_equal(got, want) + + +def test_drop_keeps_per_epoch_bounds(variable): + """Test that dropping an epoch drops its bounds too.""" + kept = variable.copy().drop([1]) + assert len(kept) == 3 + assert_allclose(kept.durations, variable.durations[[0, 2, 3]]) + + +def test_shift_time_moves_bounds_not_samples(variable): + """Test that shifting the origin does not resample anything.""" + before_lengths = [epoch.shape[1] for epoch in variable.get_data()] + before_durations = variable.durations.copy() + shifted = variable.copy().shift_time(0.1) + assert_allclose(shifted.tmin, variable.tmin + 0.1) + assert_allclose(shifted.durations, before_durations) + assert [epoch.shape[1] for epoch in shifted.get_data()] == before_lengths # -- the time axis --------------------------------------------------------- @@ -245,6 +332,30 @@ def test_fixed_epochs_still_have_times(): assert epochs.average().data.shape == (len(CH_NAMES), 71) +def test_nothing_reaches_the_user_as_an_internal_error(variable): + """Test that no public method leaks a NumPy error about lists.""" + import warnings + + names = ( + sorted(_VARIABLE_FALLBACK) + + sorted(_VARIABLE_NEEDS_POLICY) + + sorted(_VARIABLE_NOT_IMPLEMENTED) + ) + for name in names: + with warnings.catch_warnings(): + warnings.simplefilter("ignore") + try: + getattr(variable.copy(), name)() + except (NotImplementedError, RuntimeError): + pass + except TypeError as exc: + assert "argument" in str(exc), f"{name}: {exc}" + except (AttributeError, IndexError) as exc: + raise AssertionError(f"{name} leaked an internal error: {exc}") + except Exception: + pass + + # -- construction from Raw ------------------------------------------------- def _raw(n_seconds=30.0, sfreq=SFREQ): """Return a small continuous recording.""" @@ -352,3 +463,41 @@ def test_from_raw_scalar_bounds_still_scalar(): assert not epochs.variable_duration assert isinstance(epochs.tmin, float) assert epochs.get_data().shape == (2, len(CH_NAMES), 51) + + +@pytest.mark.parametrize( + "item", [slice(None, 2), slice(1, None), slice(None, None, 2), slice(None)] +) +def test_getitem_slice_selects_epochs_not_the_slice(variable, item): + """Test that slicing subsets the epochs rather than wrapping the slice.""" + want = np.arange(len(variable))[item] + subset = variable[item] + + assert len(subset) == len(want) + # a slice used to survive into the data list as a single nested element + assert all(isinstance(d, np.ndarray) and d.ndim == 2 for d in subset.get_data()) + assert_allclose(subset.durations, variable.durations[want]) + for got, idx in zip(subset.get_data(), want): + assert_array_equal(got, variable.get_data()[idx]) + + +def test_apply_function_refuses(variable): + """Test that apply_function refuses instead of indexing a list with a tuple.""" + with pytest.raises(NotImplementedError, match="not implemented"): + variable.apply_function(lambda x: x * 2) + + +def test_pick_does_not_reach_back_into_the_parent(variable): + """Test that picking replaces one object's epochs and no other's.""" + before = [epoch.shape for epoch in variable.get_data()] + + # _pick_drop_channels replaces the list contents in place, so anything + # sharing that list would be picked too + variable.copy().pick(["a", "c"]) + assert [epoch.shape for epoch in variable.get_data()] == before + + subset = variable[:2] + assert subset._data is not variable._data + subset.pick(["a"]) + assert [epoch.shape for epoch in variable.get_data()] == before + assert all(epoch.shape[0] == 1 for epoch in subset.get_data()) diff --git a/mne/utils/mixin.py b/mne/utils/mixin.py index 3addd688797..648aefa73c6 100644 --- a/mne/utils/mixin.py +++ b/mne/utils/mixin.py @@ -222,7 +222,12 @@ def _getitem( self._sanity_check_event_id() inst = self.copy() if copy else self if self._data is not None: - np.copyto(inst._data, self._data, casting="no") + if isinstance(self._data, list): + # variable-duration epochs hold one array per epoch, so there is + # nothing to copy into; copy() already produced the list + inst._data = [d.copy() for d in self._data] if copy else self._data + else: + np.copyto(inst._data, self._data, casting="no") del self select = inst._item_to_select(item) @@ -257,9 +262,21 @@ def _getitem( # will reset the index for us GetEpochsMixin.metadata.fset(inst, metadata, verbose=False) if inst.preload and select_data: - # ensure that each Epochs instance owns its own data so we can - # resize later if necessary - inst._data = np.require(inst._data[select], requirements=["O"]) + if isinstance(inst._data, list): + # `select` can still be a slice here, and iterating one yields + # the slice itself rather than the epochs it covers + inst._data = [ + inst._data[ii] for ii in np.arange(len(inst._data))[select] + ] + else: + # ensure that each Epochs instance owns its own data so we can + # resize later if necessary + inst._data = np.require(inst._data[select], requirements=["O"]) + # per-event bounds travel with the epochs they describe + if getattr(inst, "_variable_duration", False): + # an ndarray takes a slice or an index array equally well + inst._tmin_per_epoch = inst._tmin_per_epoch[select] + inst._tmax_per_epoch = inst._tmax_per_epoch[select] if drop_event_id: # update event id to reflect new content of inst inst.event_id = { @@ -767,6 +784,20 @@ def shift_time(self, tshift, relative=True) -> Self: or change the *data* values in any way. """ _check_preload(self, "shift_time") + if getattr(self, "_variable_duration", False): + # each epoch keeps its own length; only the origin moves, so this + # shifts the bounds and leaves the samples alone + if relative: + shift = tshift + else: + shift = tshift - self._tmin_per_epoch.min() + self._tmin_per_epoch = self._tmin_per_epoch + shift + self._tmax_per_epoch = self._tmax_per_epoch + shift + self._raw_times = self._raw_times + shift + self._set_times(self._raw_times) + self._update_first_last() + return self + start = tshift + (self.times[0] if relative else 0.0) new_times = start + np.arange(len(self.times)) / self.info["sfreq"] self._set_times(new_times) From 5979f12e2d52bd028aae72e04d829bcf9496de7c Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 13:57:34 +0200 Subject: [PATCH 03/11] ENH: plot variable-duration epochs without padding The epochs browser already draws its trials as a pseudo-continuous strip, concatenating them and ruling a line at each boundary, so ragged epochs need the samples they actually hold rather than a padded copy. plot() leaves the not-implemented table and becomes native. The x axis is built from those samples: lengths = per-epoch sample counts boundary_samples = np.r_[0, np.cumsum(lengths)] boundary_times = boundary_samples / sfreq n_times = boundary_samples[-1] _n_times_per_epoch returns len(times) when durations are equal, so this is one code path and reproduces the previous uniform grid exactly while never reading `times`, which variable-duration epochs refuse to provide. A window of k epochs from index i spans boundary_times[i + k] minus boundary_times[i]; _epoch_window computes that for both backends, and _get_start_stop and _load_data share _get_epoch_ix_range so the sample bounds and the concatenated array agree by construction. That is what makes the existing shape assertions meaningful for ragged windows. Arrow keys step whole epochs and shift steps whole windows, home and end ask the boundaries how many seconds an epoch is worth, and the scrollbar draws each epoch at its own width. Vertical lines mark a latency relative to each epoch's own event and are omitted from epochs too short to reach it, replacing arithmetic that took the remainder against one duration. Events map through each epoch's own window; the fixed path keeps its existing bounds, whose upper limit overshoots the last sample by |tmin|, rather than have that copied. _compute_scalings failed first of all, before any of the above, since it reshaped _data as an array. ICA sources reach the same browser without supplying the new per-epoch arrays, so those are derived from the boundaries when absent. Non-matplotlib backends decline with a message naming matplotlib until mne-qt-browser can consume boundary_times, boundary_samples and n_times, which the params dict now carries; the browser tests skip there for the same reason. --- mne/epochs.py | 1 - mne/tests/test_epochs_variable_duration.py | 7 + mne/viz/_figure.py | 101 ++++++++-- mne/viz/_mpl_figure.py | 97 ++++++--- mne/viz/epochs.py | 67 +++++-- mne/viz/tests/test_epochs.py | 219 +++++++++++++++++++++ mne/viz/utils.py | 9 +- 7 files changed, 442 insertions(+), 59 deletions(-) diff --git a/mne/epochs.py b/mne/epochs.py index a5383939fa5..94885363ea2 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -4066,7 +4066,6 @@ def _events_from_annotations(raw, events, event_id, annotations, on_missing): #: rather than a slow one, so they raise until implemented natively. _VARIABLE_NOT_IMPLEMENTED = { "filter": "filtering", - "plot": "browsing", "apply_function": "applying a function", "apply_baseline": "baseline correction", "crop": "cropping", diff --git a/mne/tests/test_epochs_variable_duration.py b/mne/tests/test_epochs_variable_duration.py index a13358f9ca9..69309e75d37 100644 --- a/mne/tests/test_epochs_variable_duration.py +++ b/mne/tests/test_epochs_variable_duration.py @@ -259,6 +259,13 @@ def test_display_methods_warn_and_fall_back(variable, meth): assert getattr(variable, meth)() is not None +def test_plot_is_not_a_fallback(): + """Test that browsing is native, not padded (see mne/viz/tests/test_epochs).""" + assert "plot" not in _VARIABLE_FALLBACK + assert "plot" not in _VARIABLE_NEEDS_POLICY + assert "plot" not in _VARIABLE_NOT_IMPLEMENTED + + # -- operations that stay native ------------------------------------------- def test_pick_keeps_durations(variable): """Test that channel selection leaves the time axis alone.""" diff --git a/mne/viz/_figure.py b/mne/viz/_figure.py index 0579820883b..1e6f324abd4 100644 --- a/mne/viz/_figure.py +++ b/mne/viz/_figure.py @@ -44,6 +44,37 @@ def __init__(self, **kwargs): vars(self).update(**kwargs) +def _epoch_window(boundary_times, start_ix, n_epochs): + """Return the start time and duration of a window of whole epochs. + + Epochs may differ in duration, so a window of ``n_epochs`` of them spans + whatever lies between the two boundaries rather than a fixed number of + seconds. ``start_ix`` is clamped so the requested epochs stay visible when + the object is long enough to allow it. + + Parameters + ---------- + boundary_times : array + Cumulative epoch edges in seconds, including both ends. + start_ix : int + Index of the first epoch to show. + n_epochs : int + Number of epochs to show. + + Returns + ------- + t_start : float + Time of the first boundary. + duration : float + Seconds spanned by the requested epochs. + """ + n_total = len(boundary_times) - 1 + n_epochs = int(np.clip(n_epochs, 1, n_total)) + start_ix = int(np.clip(start_ix, 0, n_total - n_epochs)) + stop_ix = start_ix + n_epochs + return boundary_times[start_ix], boundary_times[stop_ix] - boundary_times[start_ix] + + class BrowserBase(ABC): """A base class containing for the 2D browser. @@ -78,7 +109,13 @@ def __init__(self, **kwargs): f"Expected an instance of Raw, Epochs, or ICA, got {type(inst)}." ) - if len(inst.times) < 2: + # variable-duration epochs have no one time axis, so count the samples + # the browser will lay end to end instead + if self.mne.instance_type == "epochs": + n_inst_times = self.mne.n_times + else: + n_inst_times = len(inst.times) + if n_inst_times < 2: raise ValueError( "Data from at least two time points are required to open the browser." ) @@ -120,7 +157,12 @@ def __init__(self, **kwargs): self.mne.epoch_traces = list() self.mne.bad_epochs = list() if inst is not None: - self.mne.sampling_period = np.diff(inst.times[:2])[0] / inst.info["sfreq"] + # NB: this is 1 / sfreq**2, not a sampling period; it is only ever + # used as a small nudge before searchsorted on boundary_times, so + # keep the value while deriving it without touching inst.times + # (which variable-duration epochs refuse to provide). + sfreq = inst.info["sfreq"] + self.mne.sampling_period = (1.0 / sfreq) / sfreq # annotations self.mne.annotations = list() self.mne.hscroll_annotations = list() @@ -157,6 +199,23 @@ def __init__(self, **kwargs): self.mne.midpoints = ( np.convolve(self.mne.boundary_times, np.ones(2), mode="valid") / 2 ) + # callers that only ever deal with equal-length epochs (ICA sources) + # do not supply these, so derive them from the boundaries + sfreq = self.mne.info["sfreq"] + if not hasattr(self.mne, "boundary_samples"): + self.mne.boundary_samples = np.round( + np.asarray(self.mne.boundary_times) * sfreq + ).astype(int) + if not hasattr(self.mne, "epoch_tmins"): + n_epochs_total = len(self.mne.boundary_times) - 1 + self.mne.epoch_tmins = np.full( + n_epochs_total, float(self.mne.inst.times[0]) + ) + self.mne.epoch_tmaxs = ( + self.mne.epoch_tmins + + np.diff(self.mne.boundary_times) + - 1.0 / sfreq + ) # initialize picks and projectors self._update_picks() @@ -328,14 +387,36 @@ def _make_butterfly_selections_dict(self): # MANAGE DATA # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # + def _get_epoch_ix_range(self): + """Return the first and last+1 epoch index currently in view. + + Both :meth:`_get_start_stop` and :meth:`_load_data` go through here so + the sample bounds and the concatenated data cannot disagree, which is + what keeps the shape assertions in :meth:`_update_data` meaningful when + epochs differ in duration. + """ + # subtract one sample from tstart before searchsorted, to make sure + # we land on the left side of the boundary time (avoid precision + # errors) + ix_start = int( + np.searchsorted( + self.mne.boundary_times, self.mne.t_start - self.mne.sampling_period + ) + ) + n_total = len(self.mne.boundary_times) - 1 + ix_start = min(ix_start, max(n_total - 1, 0)) + ix_stop = min(ix_start + self.mne.n_epochs, n_total) + return ix_start, ix_stop + def _get_start_stop(self): # update time start_sec = self.mne.t_start - self.mne.first_time if self.mne.is_epochs: - start, stop = np.round( - np.array([start_sec, start_sec + self.mne.duration]) - * self.mne.info["sfreq"] - ).astype(int) + # take the samples the visible epochs really hold, so that this + # agrees with _load_data by construction rather than by arithmetic + ix_start, ix_stop = self._get_epoch_ix_range() + start = int(self.mne.boundary_samples[ix_start]) + stop = int(self.mne.boundary_samples[ix_stop]) else: # ensure our end time includes the last sample disp_duration = ( @@ -355,13 +436,7 @@ def _load_data(self, start=None, stop=None): else: return self.mne.inst[:, start:stop] else: - # subtract one sample from tstart before searchsorted, to make sure - # we land on the left side of the boundary time (avoid precision - # errors) - ix_start = np.searchsorted( - self.mne.boundary_times, self.mne.t_start - self.mne.sampling_period - ) - ix_stop = ix_start + self.mne.n_epochs + ix_start, ix_stop = self._get_epoch_ix_range() item = slice(ix_start, ix_stop) data = np.concatenate( self.mne.inst.get_data(item=item, copy=False), axis=-1 diff --git a/mne/viz/_mpl_figure.py b/mne/viz/_mpl_figure.py index ff04b3395a6..2ab30e079c9 100644 --- a/mne/viz/_mpl_figure.py +++ b/mne/viz/_mpl_figure.py @@ -58,7 +58,7 @@ from ..defaults import DEFAULTS from ..fixes import _close_event from ..utils import Bunch, _click_ch_name, logger -from ._figure import BrowserBase +from ._figure import BrowserBase, _epoch_window from .utils import ( _BLIT_KWARGS, DraggableLine, @@ -489,7 +489,7 @@ def __init__(self, inst, figsize, ica=None, xlabel="Time (s)", **kwargs): epoch_nums = self.mne.inst.selection for ix, _ in enumerate(epoch_nums): start = self.mne.boundary_times[ix] - width = np.diff(self.mne.boundary_times[:2])[0] + width = self.mne.boundary_times[ix + 1] - start ax_hscroll.add_patch( Rectangle( (start, 0), @@ -791,12 +791,20 @@ def _keypress(self, event): old_t_start = self.mne.t_start direction = 1 if key.endswith("right") else -1 if self.mne.is_epochs: - denom = 1 if key.startswith("shift") else self.mne.n_epochs + # step whole epochs, since they need not share a duration: one + # epoch normally, a whole window with shift + step = self.mne.n_epochs if key.startswith("shift") else 1 + ix_start, _ = self._get_epoch_ix_range() + self.mne.t_start, self.mne.duration = _epoch_window( + self.mne.boundary_times, + ix_start + direction * step, + self.mne.n_epochs, + ) else: denom = 1 if key.startswith("shift") else 4 - t_max = last_time - self.mne.duration - t_start = self.mne.t_start + direction * self.mne.duration / denom - self.mne.t_start = np.clip(t_start, self.mne.first_time, t_max) + t_max = last_time - self.mne.duration + t_start = self.mne.t_start + direction * self.mne.duration / denom + self.mne.t_start = np.clip(t_start, self.mne.first_time, t_max) if self.mne.t_start != old_t_start: self._update_hscroll() self._redraw(annotations=True, skip_hscroll=True) @@ -828,13 +836,17 @@ def _keypress(self, event): old_dur = self.mne.duration dur_delta = 1 if key == "end" else -1 if self.mne.is_epochs: + ix_start, _ = self._get_epoch_ix_range() # prevent from showing zero epochs, or more epochs than we have - self.mne.n_epochs = np.clip( - self.mne.n_epochs + dur_delta, 1, len(self.mne.inst) + self.mne.n_epochs = int( + np.clip(self.mne.n_epochs + dur_delta, 1, len(self.mne.inst)) + ) + # the epochs added or removed have their own durations, so ask + # the boundaries how many seconds that actually is + self.mne.t_start, new_dur = _epoch_window( + self.mne.boundary_times, ix_start, self.mne.n_epochs ) - # use the length of one epoch as duration change - min_dur = len(self.mne.inst.times) / self.mne.info["sfreq"] - new_dur = self.mne.duration + dur_delta * min_dur + min_dur = np.diff(self.mne.boundary_times).min() else: # never show fewer than 3 samples min_dur = 3 * np.diff(self.mne.inst.times[:2])[0] @@ -843,8 +855,10 @@ def _keypress(self, event): new_dur = self.mne.duration * dur_delta self.mne.duration = np.clip(new_dur, min_dur, last_time) if self.mne.duration != old_dur: - if self.mne.t_start + self.mne.duration > last_time: - self.mne.t_start = last_time - self.mne.duration + if not self.mne.is_epochs: + if self.mne.t_start + self.mne.duration > last_time: + self.mne.t_start = last_time - self.mne.duration + # (the epochs branch above already clamped t_start to a boundary) self._update_hscroll() self._redraw(annotations=True, skip_hscroll=True) elif key == "?": # help window @@ -998,7 +1012,10 @@ def _mouse_move(self, event): time = np.clip(time, self.mne.first_time, max_time) if self.mne.is_epochs: ix = np.searchsorted(self.mne.boundary_times[1:], time, side="right") - time = self.mne.boundary_times[ix] + # the epochs from here on have their own durations + time, self.mne.duration = _epoch_window( + self.mne.boundary_times, ix, self.mne.n_epochs + ) if self.mne.t_start != time: self.mne.t_start = time self._update_hscroll() @@ -1917,7 +1934,10 @@ def _check_update_hscroll_clicked(self, event): time = np.clip(time, self.mne.first_time, max_time) if self.mne.is_epochs: ix = np.searchsorted(self.mne.boundary_times[1:], time, side="right") - time = self.mne.boundary_times[ix] + # the epochs from here on have their own durations + time, self.mne.duration = _epoch_window( + self.mne.boundary_times, ix, self.mne.n_epochs + ) if self.mne.t_start != time: self.mne.t_start = time self._update_hscroll() @@ -2252,11 +2272,9 @@ def _draw_traces(self): # handle custom epoch colors (for autoreject integration) if self.mne.epoch_colors is None: # shape: n_traces × RGBA → n_traces × n_epochs × RGBA - custom_colors = np.tile( - ch_colors[:, None, :], (1, self.mne.n_epochs, 1) - ) + custom_colors = np.tile(ch_colors[:, None, :], (1, len(epoch_ix), 1)) else: - custom_colors = np.empty((len(self.mne.picks), self.mne.n_epochs, 4)) + custom_colors = np.empty((len(self.mne.picks), len(epoch_ix), 4)) for ii, _epoch_ix in enumerate(epoch_ix): this_colors = self.mne.epoch_colors[_epoch_ix] custom_colors[:, ii] = to_rgba_array( @@ -2408,26 +2426,43 @@ def _recompute_epochs_vlines(self, xdata): # special case: changed view duration w/ "home" or "end" key # (no click event, hence no xdata) if xdata is None: - xdata = np.array(self.mne.vline.get_segments())[0, 0, 0] - # compute the (continuous) times for the lines on each epoch - epoch_dur = np.diff(self.mne.boundary_times[:2])[0] - rel_time = xdata % epoch_dur - abs_time = self.mne.times[0] - xs = np.arange(self.mne.n_epochs) * epoch_dur + abs_time + rel_time - segs = np.array(self.mne.vline.get_segments()) + segments = self.mne.vline.get_segments() + if not len(segments): # no visible epoch reaches that latency + return None + xdata = np.array(segments)[0, 0, 0] + # Work out which latency relative to its own event was clicked, then + # mark that same latency on every visible epoch. Epochs need not share a + # duration, so an epoch that never reaches this latency gets no line. + sfreq = self.mne.info["sfreq"] + boundary_times = self.mne.boundary_times + clicked_ix = int( + np.clip( + np.searchsorted(boundary_times[1:], xdata, side="right"), + 0, + len(boundary_times) - 2, + ) + ) + offset = round((xdata - boundary_times[clicked_ix]) * sfreq) + latency = self.mne.epoch_tmins[clicked_ix] + offset / sfreq + ix_start, ix_stop = self._get_epoch_ix_range() + xs = list() + for ix in range(ix_start, ix_stop): + tmin, tmax = self.mne.epoch_tmins[ix], self.mne.epoch_tmaxs[ix] + if tmin - 0.5 / sfreq <= latency <= tmax + 0.5 / sfreq: + xs.append(boundary_times[ix] + (latency - tmin)) + xs = np.array(xs, float) # recreate segs from scratch in case view duration changed # (i.e., handle case when n_segments != n_epochs) segs = np.tile([[0.0], [1.0]], (len(xs), 1, 2)) # y values - segs[..., 0] = np.tile(xs[:, None], 2) # x values + segs[..., 0] = np.tile(xs[:, None], 2) if len(xs) else segs[..., 0] self.mne.vline.set_segments(segs) - return rel_time + return latency def _show_vline(self, xdata): """Show the vertical line(s).""" if self.mne.is_epochs: - # convert xdata to be epoch-relative (for the text) - rel_time = self._recompute_epochs_vlines(xdata) - xdata = rel_time + self.mne.inst.times[0] + # the label shows the latency relative to each epoch's own event + xdata = self._recompute_epochs_vlines(xdata) else: self.mne.vline.set_xdata([xdata]) self.mne.vline_hscroll.set_xdata([xdata]) diff --git a/mne/viz/epochs.py b/mne/viz/epochs.py index 669205eaf4d..3f6a4c0f90d 100644 --- a/mne/viz/epochs.py +++ b/mne/viz/epochs.py @@ -964,11 +964,24 @@ def plot_epochs( raise ValueError(msg) # handle time dimension + # + # The browser lays the epochs end to end and draws a boundary between them, + # so its x axis is "concatenated seconds". Building that axis from the + # samples each epoch really has keeps it correct when durations vary, and + # reproduces the old uniform grid exactly when they do not. n_epochs = min(n_epochs, len(epochs)) - n_times = len(epochs) * len(epochs.times) - duration = n_epochs * len(epochs.times) / sfreq + lengths = np.array( + [epochs._n_times_per_epoch(ii) for ii in range(len(epochs))], int + ) # NB: this includes start and end of data: - boundary_times = np.arange(len(epochs) + 1) * len(epochs.times) / sfreq + boundary_samples = np.concatenate([[0], np.cumsum(lengths)]) + boundary_times = boundary_samples / sfreq + n_times = int(boundary_samples[-1]) + duration = boundary_times[n_epochs] - boundary_times[0] + # one value per epoch on both paths, so the browser never needs to ask + # whether it is looking at a scalar or an array + epoch_tmins = np.broadcast_to(np.asarray(epochs.tmin, float), (len(epochs),)).copy() + epoch_tmaxs = np.broadcast_to(np.asarray(epochs.tmax, float), (len(epochs),)).copy() # events _validate_type(events, (bool, np.ndarray), "events") @@ -980,14 +993,25 @@ def plot_epochs( events = epochs.events event_nums = events[:, 2] event_samps = events[:, 0] - epoch_n_samps = len(epochs.times) + # first sample of each epoch in the raw recording it was cut from + first_samps = epochs.events[:, 0] - np.round(-epoch_tmins * sfreq).astype(int) # handle overlapping epochs (each event may show up in multiple places) - boundaries = epochs.events[:, [0]] + np.array([-1, 1]) * epochs.time_as_index( - [0, epochs.tmax] - ) - in_bounds = np.logical_and( - boundaries[:, [0]] <= event_samps, event_samps < boundaries[:, [1]] - ) + if epochs._variable_duration: + # each epoch spans the samples it actually holds; both ends are + # inclusive, matching the fixed path's effective behaviour + last_samps = first_samps + lengths - 1 + in_bounds = np.logical_and( + first_samps[:, None] <= event_samps, + event_samps <= last_samps[:, None], + ) + else: + epoch_n_samps = len(epochs.times) + boundaries = epochs.events[:, [0]] + np.array( + [-1, 1] + ) * epochs.time_as_index([0, epochs.tmax]) + in_bounds = np.logical_and( + boundaries[:, [0]] <= event_samps, event_samps < boundaries[:, [1]] + ) event_ixs = [np.nonzero(a)[0] for a in in_bounds.T] warned = False event_times = list() @@ -1000,8 +1024,13 @@ def plot_epochs( "lines may be duplicated in the plot." ) warned = True - offsets = samp - relevant_epoch_events + epochs.time_as_index(0) - this_event_times = (_ixs * epoch_n_samps + offsets) / sfreq + if epochs._variable_duration: + # position inside each containing epoch, then along the strip + offsets = samp - first_samps[_ixs] + this_event_times = boundary_times[_ixs] + offsets / sfreq + else: + offsets = samp - relevant_epoch_events + epochs.time_as_index(0) + this_event_times = (_ixs * epoch_n_samps + offsets) / sfreq event_times.extend(this_event_times) event_numbers.extend([num] * len(_ixs)) event_nums = np.array(event_numbers) @@ -1059,6 +1088,9 @@ def plot_epochs( time_format="float", decim=decim, boundary_times=boundary_times, + boundary_samples=boundary_samples, + epoch_tmins=epoch_tmins, + epoch_tmaxs=epoch_tmaxs, # events event_id_rev=event_id_rev, event_color_dict=event_color_dict, @@ -1100,6 +1132,17 @@ def plot_epochs( figure_class=figure_class, ) + if epochs._variable_duration: + from ._figure import get_browser_backend + + backend_name = get_browser_backend() + if backend_name != "matplotlib": + raise NotImplementedError( + f"Browsing variable-duration epochs is not implemented for the " + f"{backend_name} backend yet, only for matplotlib. Select it " + 'with mne.viz.set_browser_backend("matplotlib").' + ) + fig = _get_browser(show=show, block=block, **params) return fig diff --git a/mne/viz/tests/test_epochs.py b/mne/viz/tests/test_epochs.py index 874efa88e77..78fbcf04513 100644 --- a/mne/viz/tests/test_epochs.py +++ b/mne/viz/tests/test_epochs.py @@ -3,10 +3,12 @@ # Copyright the MNE-Python contributors. import platform +import warnings import matplotlib.pyplot as plt import numpy as np import pytest +from numpy.testing import assert_allclose, assert_array_equal from mne import Epochs, EpochsArray, create_info from mne.datasets import testing @@ -511,3 +513,220 @@ def test_plot_epochs_selection_butterfly(raw, browser_backend): epochs = Epochs(raw, events, tmin=0, tmax=0.5, preload=True, baseline=None) assert len(epochs) == 1 epochs.plot(group_by="selection", butterfly=True) + + +# -- variable-duration epochs ---------------------------------------------- +SFREQ_VAR = 100.0 +LENGTHS_VAR = (100, 250, 75, 180) # deliberately very different + + +def _variable_epochs(tmin=None): + """Return epochs whose trials have deliberately unequal lengths.""" + n = len(LENGTHS_VAR) + info = create_info(["a", "b", "c"], SFREQ_VAR, "eeg") + rng = np.random.default_rng(0) + data = [rng.standard_normal((3, length)) * 1e-6 for length in LENGTHS_VAR] + events = np.column_stack( + [np.arange(n) * 1000 + 500, np.zeros(n, int), np.ones(n, int)] + ) + tmin = np.zeros(n) if tmin is None else np.asarray(tmin, float) + return EpochsArray( + data, + info, + events=events, + tmin=tmin, + event_id={"x": 1}, + baseline=None, + verbose=False, + ) + + +def _boundaries(): + """Return the boundary times the browser should derive.""" + return np.concatenate([[0], np.cumsum(LENGTHS_VAR)]) / SFREQ_VAR + + +@pytest.fixture +def variable_epochs(): + """Epochs of unequal duration sharing tmin=0.""" + return _variable_epochs() + + +def test_plot_variable_duration_is_native(variable_epochs, browser_backend): + """Test that browsing ragged epochs neither warns nor pads.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + + def _boom(*args, **kwargs): + raise AssertionError("plot() fell back to as_fixed() instead of browsing") + + variable_epochs.as_fixed = _boom + with warnings.catch_warnings(): + warnings.simplefilter("error") # any warning fails the test + fig = variable_epochs.plot(n_epochs=2) + assert not np.isnan(fig.mne.data).any() + + +def test_plot_variable_duration_boundaries(variable_epochs, browser_backend): + """Test that the browser lays epochs end to end at their true lengths.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + fig = variable_epochs.plot(n_epochs=2) + assert_allclose(fig.mne.boundary_times, _boundaries()) + assert_array_equal(fig.mne.boundary_samples, np.r_[0, np.cumsum(LENGTHS_VAR)]) + # the concatenated axis holds every sample once, with nothing invented + assert fig.mne.n_times == sum(LENGTHS_VAR) + assert fig.mne.n_times != len(LENGTHS_VAR) * max(LENGTHS_VAR) + + +def test_plot_variable_duration_window_spans_whole_epochs( + variable_epochs, browser_backend +): + """Test that n_epochs means epochs, not a representative duration.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + boundaries = _boundaries() + for n_epochs in (1, 2, 3, 4): + fig = variable_epochs.plot(n_epochs=n_epochs) + assert fig.mne.duration == pytest.approx(boundaries[n_epochs]) + assert fig.mne.data.shape[-1] == sum(LENGTHS_VAR[:n_epochs]) + + +def test_plot_variable_duration_data_is_unpadded(variable_epochs, browser_backend): + """Test that a view holds exactly the source samples, in order.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + fig = variable_epochs.plot(n_epochs=2) + source = variable_epochs.get_data() + for keys in ([], ["right"], ["right", "right"]): + for key in keys: + fig._fake_keypress(key) + start, stop = fig._get_start_stop() + data, times = fig._load_data(start, stop) + ix_start, ix_stop = fig._get_epoch_ix_range() + want = np.concatenate(source[ix_start:ix_stop], axis=-1) + # the raw window is the source samples, in order, with nothing added + assert_array_equal(data, want) + assert not np.isnan(data).any() + # and the sample bounds agree with it, which is what lets the shape + # assertions in _update_data mean something for ragged windows + assert data.shape[-1] == stop - start + assert len(times) == stop - start + + +def test_plot_variable_duration_navigation(variable_epochs, browser_backend): + """Test that arrow keys move by epochs and land on real boundaries.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + boundaries = _boundaries() + fig = variable_epochs.plot(n_epochs=2) + assert fig.mne.t_start == pytest.approx(boundaries[0]) + + fig._fake_keypress("right") # one epoch + assert fig.mne.t_start == pytest.approx(boundaries[1]) + assert fig.mne.duration == pytest.approx(boundaries[3] - boundaries[1]) + + fig._fake_keypress("right") # clamped: two epochs must stay visible + assert fig.mne.t_start == pytest.approx(boundaries[2]) + assert fig.mne.duration == pytest.approx(boundaries[4] - boundaries[2]) + + fig._fake_keypress("left") + assert fig.mne.t_start == pytest.approx(boundaries[1]) + + # shift moves a whole window + fig._fake_keypress("shift+left") + assert fig.mne.t_start == pytest.approx(boundaries[0]) + fig._fake_keypress("shift+right") + assert fig.mne.t_start == pytest.approx(boundaries[2]) + + +def test_plot_variable_duration_home_end(variable_epochs, browser_backend): + """Test that home/end change the epoch count and recompute the duration.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + boundaries = _boundaries() + fig = variable_epochs.plot(n_epochs=2) + assert fig.mne.duration == pytest.approx(boundaries[2]) + + fig._fake_keypress("end") # show one more epoch + assert fig.mne.n_epochs == 3 + assert fig.mne.duration == pytest.approx(boundaries[3]) + # the epoch added is 75 samples, not a repeat of the first one + assert fig.mne.duration != pytest.approx(boundaries[2] * 3 / 2) + + fig._fake_keypress("home") + assert fig.mne.n_epochs == 2 + assert fig.mne.duration == pytest.approx(boundaries[2]) + + +def test_plot_variable_duration_hscroll_patches(variable_epochs, browser_backend): + """Test that the scrollbar draws each epoch at its own width.""" + if browser_backend.name != "matplotlib": + pytest.skip("scrollbar patches are matplotlib-specific") + fig = variable_epochs.plot(n_epochs=2) + widths = [p.get_width() for p in fig.mne.ax_hscroll.patches[: len(LENGTHS_VAR)]] + assert_allclose(widths, np.diff(_boundaries())) + assert len(set(np.round(widths, 6))) == len(LENGTHS_VAR) # all different + + +def test_plot_variable_duration_bad_epoch(variable_epochs, browser_backend): + """Test that a click finds the right epoch when the widths differ.""" + if browser_backend.name != "matplotlib": + pytest.skip("variable-duration browsing is matplotlib-only") + boundaries = _boundaries() + fig = variable_epochs.plot(n_epochs=4) + y = fig.mne.traces[0].get_ydata()[0] + # click inside epoch 2, which starts well past twice the first epoch's width + x = (boundaries[2] + boundaries[3]) / 2 + fig._fake_click((x, y), xform="data") + assert list(fig.mne.bad_epochs) == [variable_epochs.selection[2]] + fig._fake_click((x, y), xform="data") # unmark + assert list(fig.mne.bad_epochs) == [] + + +def test_plot_variable_duration_vline_latency(variable_epochs, browser_backend): + """Test that vlines mark a latency, skipping epochs that never reach it.""" + if browser_backend.name != "matplotlib": + pytest.skip("vline segments are matplotlib-specific") + boundaries = _boundaries() + fig = variable_epochs.plot(n_epochs=4) + + # 1.5 s exists only in the 250- and 180-sample epochs + latency = 1.5 + fig._fake_click((boundaries[1] + latency, 0.5), xform="data") + xs = np.sort(np.array(fig.mne.vline.get_segments())[:, 0, 0]) + assert_allclose(xs, [boundaries[1] + latency, boundaries[3] + latency]) + + # 0.5 s exists in every epoch + fig._fake_click((boundaries[0] + 0.5, 0.5), xform="data") + xs = np.sort(np.array(fig.mne.vline.get_segments())[:, 0, 0]) + assert_allclose(xs, boundaries[:4] + 0.5) + + +def test_plot_variable_duration_events(browser_backend): + """Test that events map into the strip using each epoch's own window.""" + if browser_backend.name != "matplotlib": + pytest.skip("event line segments are matplotlib-specific") + # unequal tmin as well as unequal duration + tmin = np.array([0.0, -0.2, 0.0, -0.5]) + epochs = _variable_epochs(tmin=tmin) + boundaries = _boundaries() + first_samps = epochs.events[:, 0] - np.round(-tmin * SFREQ_VAR).astype(int) + + fig = epochs.plot(n_epochs=4, events=True) + lines, _ = _get_event_lines_and_texts(fig) + got = np.sort(np.array(lines)[:, 0, 0]) + # each defining event sits at its own offset inside its own epoch + want = np.sort( + boundaries[: len(LENGTHS_VAR)] + (epochs.events[:, 0] - first_samps) / SFREQ_VAR + ) + assert_allclose(got, want) + + +def test_plot_variable_duration_refuses_other_backends(variable_epochs, monkeypatch): + """Test that non-matplotlib backends decline rather than fail obscurely.""" + import mne.viz._figure + + monkeypatch.setattr(mne.viz._figure, "get_browser_backend", lambda: "qt") + with pytest.raises(NotImplementedError, match="not implemented for the qt"): + variable_epochs.plot() diff --git a/mne/viz/utils.py b/mne/viz/utils.py index 80a5b146ec8..b36440fe0f9 100644 --- a/mne/viz/utils.py +++ b/mne/viz/utils.py @@ -1401,7 +1401,8 @@ def _compute_scalings(scalings, inst, remove_dc=False, duration=10): data = inst._read_segment(smin, smax) elif isinstance(inst, BaseEpochs): # Load a random subset of epochs up to 100mb in size - n_epochs = 1e8 // (len(inst.ch_names) * len(inst.times) * 8) + longest = max(inst._n_times_per_epoch(ii) for ii in range(len(inst))) + n_epochs = 1e8 // (len(inst.ch_names) * longest * 8) n_epochs = int(np.clip(n_epochs, 1, len(inst))) ixs_epochs = np.random.default_rng(0).choice( len(inst), n_epochs, replace=False @@ -1410,7 +1411,11 @@ def _compute_scalings(scalings, inst, remove_dc=False, duration=10): else: data = inst._data if isinstance(inst, BaseEpochs): - data = inst._data.swapaxes(0, 1).reshape([len(inst.ch_names), -1]) + if isinstance(inst._data, list): + # variable-duration epochs: one array per epoch, already channels-first + data = np.concatenate(inst._data, axis=-1) + else: + data = inst._data.swapaxes(0, 1).reshape([len(inst.ch_names), -1]) # Iterate through ch types and update scaling if ' auto' for key, value in scalings.items(): if key not in ch_types or value != "auto": From bed56d5b7b13b631479ab67d3ed3118ce3cd8f5d Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 13:57:48 +0200 Subject: [PATCH 04/11] DOC: demonstrate variable-duration epochs Builds epochs straight from the Sleep Physionet hypnogram durations, so no new dataset is needed and tools/circleci_download.sh already prefetches it. Bouts over five minutes are set aside to keep the padded array small, leaving 130 epochs from 30 to 300 s. Walks through what the object holds, which operations are unaffected by ragged trials, which ones refuse and why, browsing them at their own lengths, and what as_fixed() reports: 130 epochs at t=0 falling to 1 by 300 s, which is the reason average() cannot return an ordinary Evoked. The browsing section picks bouts by taking the first occurrence of each distinct value in `durations`, so five different lengths are guaranteed rather than hoped for; here that is 120, 30, 150, 60 and 90 s. It runs under use_browser_backend("matplotlib"), as several other tutorials already do, because the doc build exports MNE_BROWSER_BACKEND=qt and qt is tried first, and the PyQtGraph backend does not handle ragged epochs yet. The closing section points at the sleep-staging tutorial, where fixed 30 s windows are the right representation, so the two are not read as alternatives. --- doc/changes/dev/14210.newfeature.rst | 2 +- .../epochs/70_variable_duration_epochs.py | 304 ++++++++++++++++++ 2 files changed, 305 insertions(+), 1 deletion(-) create mode 100644 tutorials/epochs/70_variable_duration_epochs.py diff --git a/doc/changes/dev/14210.newfeature.rst b/doc/changes/dev/14210.newfeature.rst index fc0cde80c82..a657833cd9c 100644 --- a/doc/changes/dev/14210.newfeature.rst +++ b/doc/changes/dev/14210.newfeature.rst @@ -1 +1 @@ -Allow :class:`mne.Epochs` to hold trials of different duration by passing ``tmin`` and/or ``tmax`` as arrays with one entry per event, with :meth:`mne.Epochs.as_fixed` to obtain a fixed-duration copy spanning their union together with the number of epochs contributing at each time point, by `Sina Esmaeili`_. +Allow :class:`mne.Epochs` to hold trials of different duration by passing ``tmin`` and/or ``tmax`` as arrays with one entry per event, with :meth:`mne.Epochs.as_fixed` to obtain a fixed-duration copy spanning their union together with the number of epochs contributing at each time point, and with :meth:`mne.Epochs.plot` browsing them at their own lengths rather than padding, demonstrated in :ref:`tut-variable-duration-epochs`, by `Sina Esmaeili`_. diff --git a/tutorials/epochs/70_variable_duration_epochs.py b/tutorials/epochs/70_variable_duration_epochs.py new file mode 100644 index 00000000000..52afe8bd2e5 --- /dev/null +++ b/tutorials/epochs/70_variable_duration_epochs.py @@ -0,0 +1,304 @@ +""" +.. _tut-variable-duration-epochs: + +Epochs whose trials have different durations +============================================ + +Most epoching starts from an event and takes the same window around every one of +them, which gives a rectangular ``(n_epochs, n_channels, n_times)`` array and one +time axis shared by every trial. Some experiments do not fit that shape. A gait +cycle, a spoken word, a reaching movement and a sleep stage all last as long as +they last, and the duration is often the thing being studied. + +The usual way to handle this is to pick a fixed window and accept the +consequences: a window long enough for the longest trial pads the short ones, and +a window short enough for the shortest one truncates the rest. This tutorial +shows the other option, keeping each trial at the length it actually had, and +what the resulting object can and cannot do. + +We use the Sleep Physionet data, where the hypnogram annotations mark sleep stage +bouts and each bout carries its own duration. +""" + +# Authors: The MNE-Python contributors. +# License: BSD-3-Clause +# Copyright the MNE-Python contributors. + +# %% + +import matplotlib.pyplot as plt +import numpy as np + +import mne +from mne.datasets.sleep_physionet.age import fetch_data + +psg_file, hypnogram_file = fetch_data(subjects=[0], recording=[1])[0] + +raw = mne.io.read_raw_edf( + psg_file, + stim_channel=False, + preload=True, + verbose="error", # ignore issues with stored filter settings +) +raw.pick(["EEG Fpz-Cz", "EEG Pz-Oz"]) + +annotations = mne.read_annotations(hypnogram_file) +raw.set_annotations(annotations, emit_warning=False) + +# %% +# The annotations already carry durations +# --------------------------------------- +# +# Each hypnogram entry marks one bout of a sleep stage, and +# :class:`~mne.Annotations` stores its ``duration`` alongside its ``onset``. That +# duration is not a constant. + +stages = { + "Sleep stage 1": 1, + "Sleep stage 2": 2, + "Sleep stage 3": 3, + "Sleep stage 4": 3, # stages 3 and 4 are conventionally merged + "Sleep stage R": 4, +} +event_id = {"N1": 1, "N2": 2, "N3/4": 3, "REM": 4} + +onset = annotations.onset +duration = annotations.duration +description = np.array(annotations.description) + +keep = np.array([desc in stages for desc in description]) +# a five minute cap keeps the padded array in the last section small; nothing +# about the container requires it +keep &= duration <= 300.0 + +onset, duration = onset[keep], duration[keep] +description = description[keep] + +print(f"{len(duration)} bouts, {duration.min():.0f} to {duration.max():.0f} s") +print(f"median {np.median(duration):.0f} s") + +# %% +# Building epochs that keep those durations +# ----------------------------------------- +# +# ``tmin`` and ``tmax`` accept one value per event as well as a single number. +# Here every bout starts at its own onset, so ``tmin`` is zero throughout and +# ``tmax`` is the bout's own length. The last sample is included, which is why +# ``tmax`` is one sample short of the full duration. + +sfreq = raw.info["sfreq"] +events = np.column_stack( + [ + np.round((onset - raw.first_time) * sfreq).astype(int), + np.zeros(len(onset), int), + np.array([stages[desc] for desc in description]), + ] +) + +epochs = mne.Epochs( + raw, + events, + event_id, + tmin=np.zeros(len(events)), + tmax=duration - 1.0 / sfreq, + baseline=None, + preload=True, +) +print(epochs) + +# %% +# The object reports that its trials are not all the same length, and the +# durations it holds are the ones the annotations described. + +print(f"variable_duration: {epochs.variable_duration}") +print(f"durations: {epochs.durations.min():.0f} to {epochs.durations.max():.0f} s") + +# %% +# Because bounds that carry no variation collapse back to a single value, this +# only changes behaviour when the durations really do differ. Passing equal +# bounds gives an ordinary fixed-duration ``Epochs``. + +n2_events = events[events[:, 2] == event_id["N2"]][:5] +fixed_bounds = mne.Epochs( + raw, + n2_events, + {"N2": event_id["N2"]}, + tmin=np.zeros(len(n2_events)), + tmax=np.full(len(n2_events), 29.99), + baseline=None, + preload=True, + verbose=False, +) +print(f"equal bounds -> variable_duration: {fixed_bounds.variable_duration}") + +# %% +# Getting the data out +# -------------------- +# +# There is no rectangular array to return, so :meth:`~mne.Epochs.get_data` gives +# a list with one ``(n_channels, n_times)`` array per epoch. Nothing is padded +# and nothing is cut: each array holds exactly the samples that the bout covered +# in the continuous recording. + +data = epochs.get_data() +print(f"{len(data)} arrays, first four shapes {[d.shape for d in data[:4]]}") + +lengths = np.array([d.shape[-1] for d in data]) +print(f"total samples held: {lengths.sum()}") +print(f"a rectangular array would hold: {lengths.max() * len(lengths)}") + +# %% +# For the same reason there is no single ``times`` attribute. Each epoch has its +# own time axis, which :meth:`~mne.Epochs.get_times` returns. + +for idx in (0, 1): + t = epochs.get_times(idx) + print(f"epoch {idx}: {len(t)} samples, {t[0]:.2f} to {t[-1]:.2f} s") + +# %% +# Asking for ``epochs.times`` raises rather than inventing an axis. Returning the +# longest epoch's axis would make ``len(epochs.times)`` disagree with the data +# for every other epoch while looking perfectly normal. + +try: + epochs.times +except RuntimeError as err: + print(f"RuntimeError: {err}") + +# %% +# Operations that do not touch the time axis +# ------------------------------------------ +# +# Selecting epochs, selecting channels and dropping epochs all work as usual, +# because none of them care how long each trial is. The per-epoch bounds travel +# with the epochs they belong to. + +n2 = epochs["N2"] +print( + f"epochs['N2']: {len(n2)} epochs, " + f"{n2.durations.min():.0f} to {n2.durations.max():.0f} s" +) + +first_ten = epochs[:10] +print(f"epochs[:10]: durations {first_ten.durations.round(0)}") + +one_channel = epochs.copy().pick(["EEG Pz-Oz"]) +print( + f"after pick: {one_channel.ch_names}, durations unchanged: " + f"{np.array_equal(one_channel.durations, epochs.durations)}" +) + +# %% +# Browsing them +# ------------- +# +# :meth:`~mne.Epochs.plot` shows each bout at the length it really has. The +# browser lays the variable-length blocks end to end and rules a line between +# them, so the vertical boundaries are unevenly spaced: a 30 second bout takes a +# fifth of the width of a 150 second one. Nothing is padded or truncated to make +# the picture rectangular, and :meth:`~mne.Epochs.as_fixed` is not involved. +# +# Pick a handful of bouts with genuinely different lengths, taking the first +# occurrence of each distinct duration rather than trusting the first few epochs +# to differ. + +_, first_of_each = np.unique(epochs.durations, return_index=True) +browse_idx = np.sort(first_of_each[:5]) +browse_epochs = epochs[browse_idx] +print(f"browsing durations: {browse_epochs.durations.round(0)} s") + +# the browser's time axis is the real samples, laid end to end +n_browser_samples = sum( + len(browse_epochs.get_times(ii)) for ii in range(len(browse_epochs)) +) +print(f"{n_browser_samples} samples in total, none of them padding") + +# %% +# Browsing variable-duration epochs currently needs the Matplotlib backend; the +# PyQtGraph one does not handle ragged epochs yet. + +with mne.viz.use_browser_backend("matplotlib"): + browse_epochs.plot(n_epochs=len(browse_epochs), picks="eeg") + +# %% +# Operations that need one time axis +# ---------------------------------- +# +# Averaging is the clearest case. :class:`~mne.Evoked` holds one array and one +# ``nave``, and there is no honest way to fill either when the trials stop at +# different times. Rather than pad quietly, the reduction refuses and says what +# it would need. + +try: + epochs.average() +except NotImplementedError as err: + print(f"NotImplementedError: {err}") + +# %% +# Making the padding explicit +# --------------------------- +# +# When a rectangular array is genuinely what you want, +# :meth:`~mne.Epochs.as_fixed` produces one. It returns the padded +# :class:`~mne.EpochsArray` together with the number of epochs contributing at +# each sample, so the cost of the padding is visible rather than implied. + +padded, n_contributing = epochs.as_fixed() +print(f"padded shape: {padded.get_data().shape}") +print( + f"contributing: {n_contributing.max()} at the start, " + f"{n_contributing.min()} at the end" +) + +held = lengths.sum() * len(epochs.ch_names) +allocated = padded.get_data().size +print(f"padding waste: {100 * (1 - held / allocated):.1f}%") + +# %% +# That second return value is the point of the method. Plotted against time it +# shows how quickly the epochs stop contributing, which is exactly the +# information an averaged :class:`~mne.Evoked` cannot carry. + +fig, ax = plt.subplots(figsize=(8, 4), layout="constrained") +times = padded.times +ax.fill_between(times, n_contributing, step="post", alpha=0.25) +ax.plot(times, n_contributing, drawstyle="steps-post") + +half = len(epochs) / 2 +crossing = times[np.argmax(n_contributing < half)] +ax.axhline(half, color="0.4", ls=":", lw=1) +ax.axvline(crossing, color="0.4", ls=":", lw=1) +ax.annotate( + f"half the epochs have ended by {crossing:.0f} s", + xy=(crossing, half), + xytext=(crossing + 20, len(epochs) * 0.7), + arrowprops=dict(arrowstyle="->", color="0.4"), +) + +ax.set( + xlabel="Time (s)", + ylabel="Epochs contributing", + title="How many sleep-stage bouts are still running", + xlim=(0, times[-1]), + ylim=(0, len(epochs) * 1.05), +) + +# %% +# Reading the figure from left to right: every bout contributes at the start, +# and by the end a single long bout is holding up the whole window. An average +# over this padded array would combine all of them at ``t = 0`` and one of them +# at the right-hand edge, while reporting one ``nave`` for the lot. Keeping the +# count alongside the data is what makes that visible. +# +# When fixed windows are the right choice +# --------------------------------------- +# +# None of this argues against fixed-length epochs. Sleep staging is a good +# example of when they are correct: :ref:`tut-sleep-stage-classif` classifies 30 +# second windows, so it passes ``chunk_duration=30.`` to +# :func:`mne.events_from_annotations` and deliberately turns each bout into a +# series of equal windows. That is the right representation when the window is +# the unit of analysis. +# +# Variable-duration epochs are for the other case, when the bout itself is the +# unit and its length is part of what is being measured. From 0cffb89615650cf194031b53366898a287a989e3 Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 15:33:48 +0200 Subject: [PATCH 05/11] ENH: crop variable-duration epochs Cropping asks for a window in seconds, and that question has an answer for each trial on its own: keep the samples inside it. No epoch has to be padded, stretched or compared with any other, so crop leaves the not-implemented table. The requested window is applied to every epoch independently and clamped to that epoch's own bounds where it reaches past them, which is what the fixed path does against its single interval. Selections for every epoch are computed before anything is written, so a window that misses one epoch fails and leaves the object as it was rather than dropping it. That failure comes from _time_mask seeing an inverted interval once tmax has been clamped back, which is the same route the fixed path takes. Clamping is reported once per bound rather than once per epoch, and only when it happened. A clamped tmax keeps that epoch's last sample even when include_tmax is False, matching the fixed path. Bounds are taken from the samples that survived, never from the requested float, so len(get_times(i)) continues to describe the block. Cropping can also remove the variation: when every epoch ends up on the same axis the blocks are stacked and the object becomes an ordinary Epochs again, which is checked by sample index and length rather than by comparing floats. The reductions then return on their own, since the wrappers ask about _variable_duration when they are called. ExtendedTimeMixin is untouched. It is shared with Raw, Evoked and TFR, and this behaviour belongs to Epochs. --- mne/epochs.py | 115 +++++++++++++- mne/tests/test_epochs_variable_duration.py | 176 +++++++++++++++++++++ 2 files changed, 290 insertions(+), 1 deletion(-) diff --git a/mne/epochs.py b/mne/epochs.py index 94885363ea2..556dc8c8f1d 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -97,6 +97,7 @@ _prepare_read_metadata, _prepare_write_metadata, _scale_dataframe_data, + _time_mask, _validate_type, check_fname, check_random_state, @@ -2729,6 +2730,10 @@ def crop( # XXX this could be made to work on non-preloaded data... _check_preload(self, "Modifying data of epochs") + if self._variable_duration: + self._crop_variable(tmin, tmax, include_tmax) + return self + super().crop(tmin=tmin, tmax=tmax, include_tmax=include_tmax) # Adjust rejection period @@ -2746,6 +2751,115 @@ def crop( self.reject_tmax = self.tmax return self + def _crop_variable(self, tmin, tmax, include_tmax): + """Crop each epoch on its own time axis. + + The requested window is a physical interval in seconds, so it is applied + to every epoch independently and clamped to the epoch's own bounds when + it reaches past them. Nothing is padded, interpolated or aligned; an + epoch that the window misses entirely makes the whole call fail, the + same way it would for that epoch on its own. + + Parameters + ---------- + tmin : float | None + Start of the window, or ``None`` for each epoch's own start. + tmax : float | None + End of the window, or ``None`` for each epoch's own end. + include_tmax : bool + Whether to keep the sample at ``tmax``. + """ + for name in ("reject_tmin", "reject_tmax"): + if getattr(self, name, None) is not None: + raise NotImplementedError( + f"{name} is not implemented for variable-duration epochs, " + "because the window is not guaranteed to exist in every " + "epoch." + ) + + sfreq = float(self.info["sfreq"]) + # First pass: work out every selection while changing nothing, so that a + # window that misses one epoch leaves the object as it was. + masks = list() + clamped_tmin = clamped_tmax = False + for ii in range(len(self.events)): + times = self.get_times(ii) + this_tmin, this_tmax = tmin, tmax + this_include_tmax = include_tmax + if this_tmin is None: + this_tmin = times[0] + elif this_tmin < times[0]: + clamped_tmin = True + this_tmin = times[0] + if this_tmax is None: + this_tmax = times[-1] + elif this_tmax > times[-1]: + clamped_tmax = True + this_tmax = times[-1] + # matches the fixed path: a clamped end keeps its last sample + this_include_tmax = True + # _time_mask raises when the window is inverted, which is what an + # entirely out-of-range request collapses to once tmax is clamped + mask = _time_mask( + times, + this_tmin, + this_tmax, + sfreq=sfreq, + include_tmax=this_include_tmax, + ) + if not mask.any(): + raise ValueError( + f"tmin ({tmin}) and tmax ({tmax}) would leave epoch {ii} " + f"with no samples; it spans {times[0]:g} to {times[-1]:g} s." + ) + masks.append(mask) + + # One warning per bound however many epochs needed clamping + if clamped_tmin: + warn( + "tmin is not in time interval for every epoch. tmin is set to " + "each of those epochs' own first sample." + ) + if clamped_tmax: + warn( + "tmax is not in time interval for every epoch. tmax is set to " + "each of those epochs' own last sample." + ) + + # Second pass: apply + ragged = self._data + assert ragged is not None # variable-duration epochs are always preloaded + starts = np.empty(len(masks)) + stops = np.empty(len(masks)) + for ii, mask in enumerate(masks): + kept = self.get_times(ii)[mask] + ragged[ii] = ragged[ii][..., mask] + starts[ii], stops[ii] = kept[0], kept[-1] + self._tmin_per_epoch = starts + self._tmax_per_epoch = stops + + # Cropping can remove the variation entirely. Compare the axes by their + # sample index and length rather than by float equality: every epoch is + # regularly sampled at one sfreq, so that pair identifies an axis + # exactly and does not depend on how the bound was rounded. + first_idx = np.round(starts * sfreq).astype(int) + lengths = np.array([epoch.shape[-1] for epoch in ragged]) + if ( + len(masks) + and (first_idx == first_idx[0]).all() + and (lengths == lengths[0]).all() + ): + self._data = np.stack(list(ragged)) + self._variable_duration = False + self._tmin_per_epoch = None # ty: ignore[invalid-assignment] + self._tmax_per_epoch = None # ty: ignore[invalid-assignment] + start_idx, stop_idx = first_idx[0], first_idx[0] + lengths[0] - 1 + else: + start_idx = int(round(starts.min() * sfreq)) + stop_idx = int(round(stops.max() * sfreq)) + self._raw_times = np.arange(start_idx, stop_idx + 1) / sfreq + self._set_times(self._raw_times) + def copy(self) -> Self: """Return copy of Epochs instance. @@ -4068,7 +4182,6 @@ def _events_from_annotations(raw, events, event_id, annotations, on_missing): "filter": "filtering", "apply_function": "applying a function", "apply_baseline": "baseline correction", - "crop": "cropping", "decimate": "decimation", "resample": "resampling", "save": "writing to FIF", diff --git a/mne/tests/test_epochs_variable_duration.py b/mne/tests/test_epochs_variable_duration.py index 69309e75d37..b38c8b0e3ea 100644 --- a/mne/tests/test_epochs_variable_duration.py +++ b/mne/tests/test_epochs_variable_duration.py @@ -4,6 +4,8 @@ # License: BSD-3-Clause # Copyright the MNE-Python contributors. +import warnings + import numpy as np import pytest from numpy.testing import assert_allclose, assert_array_equal @@ -508,3 +510,177 @@ def test_pick_does_not_reach_back_into_the_parent(variable): subset.pick(["a"]) assert [epoch.shape for epoch in variable.get_data()] == before assert all(epoch.shape[0] == 1 for epoch in subset.get_data()) + + +# -- crop ------------------------------------------------------------------- +def _crop_oracle(epochs, idx, **kwargs): + """Crop epoch ``idx`` as an ordinary one-epoch Epochs, for comparison.""" + data = epochs.get_data()[idx] + tmin = np.atleast_1d(epochs.tmin)[idx] + one = EpochsArray( + data[None], + create_info(list(epochs.ch_names), SFREQ, "eeg"), + tmin=float(tmin), + baseline=None, + verbose=False, + ) + return one.crop(**kwargs) + + +def _assert_crop_matches_mne(epochs, **kwargs): + """Assert every cropped epoch equals ordinary MNE cropping it alone.""" + wanted = [_crop_oracle(epochs, ii, **kwargs) for ii in range(len(epochs))] + got = epochs.copy().crop(**kwargs) + data = got.get_data() + for ii, one in enumerate(wanted): + assert_allclose(data[ii], one.get_data()[0]) + times = got.get_times(ii) if got.variable_duration else got.times + assert_allclose(times, one.times) + return got + + +@pytest.mark.parametrize( + "kwargs", + [ + dict(tmin=0.0, tmax=0.4), # both bounds + dict(tmin=0.1), # only tmin + dict(tmax=0.45), # only tmax + dict(tmin=-0.05, tmax=0.45, include_tmax=True), + dict(tmin=-0.05, tmax=0.45, include_tmax=False), + ], +) +def test_crop_matches_mne_per_epoch(variable, kwargs): + """Test that cropping equals ordinary MNE applied to each epoch alone.""" + _assert_crop_matches_mne(variable, **kwargs) + + +def test_crop_matches_mne_with_unequal_tmin(kwargs=None): + """Test parity when both bounds differ between epochs.""" + epochs = _make([-0.2, -0.35, 0.0, -0.1], [0.5, 0.9, 0.7, 0.6]) + _assert_crop_matches_mne(epochs, tmin=0.05, tmax=0.4) + _assert_crop_matches_mne(epochs, tmax=0.5) + + +def test_crop_keeps_the_object_ragged(variable): + """Test that unequal durations survive a crop that does not equalise them.""" + cropped = variable.copy().crop(tmin=0.0) + assert cropped.variable_duration + assert len(np.unique(cropped.durations)) > 1 + assert isinstance(cropped._data, list) + + +def test_crop_clamps_each_epoch_and_warns_once(variable): + """Test that a bound past some epochs clamps per epoch, warning once.""" + before = variable.durations.copy() + with pytest.warns(RuntimeWarning, match="tmax is not in time interval") as rec: + cropped = variable.copy().crop(tmax=99.0) + assert len(rec) == 1 # not one per epoch + # each epoch kept everything it had + assert_allclose(cropped.durations, before) + assert cropped.variable_duration + + with pytest.warns(RuntimeWarning, match="tmin is not in time interval") as rec: + cropped = variable.copy().crop(tmin=-99.0) + assert len(rec) == 1 + assert_allclose(cropped.durations, before) + + +def test_crop_does_not_warn_when_nothing_is_clamped(variable): + """Test that a window inside every epoch is silent.""" + with warnings.catch_warnings(): + warnings.simplefilter("error") + variable.copy().crop(tmin=0.0, tmax=0.4) + + +def test_crop_clamped_tmax_keeps_the_last_sample(variable): + """Test that clamping tmax includes that epoch's final sample.""" + lengths = [epoch.shape[-1] for epoch in variable.get_data()] + with pytest.warns(RuntimeWarning, match="tmax is not in time interval"): + # include_tmax=False must not drop the endpoint that clamping produced + cropped = variable.copy().crop(tmax=99.0, include_tmax=False) + assert [epoch.shape[-1] for epoch in cropped.get_data()] == lengths + + +def test_crop_outside_every_sample_fails_cleanly(variable): + """Test that a window missing an epoch refuses and changes nothing.""" + before_data = [epoch.copy() for epoch in variable.get_data()] + before_tmin = np.array(variable.tmin) + before_tmax = np.array(variable.tmax) + with pytest.raises(ValueError, match="must be less than or equal to"): + variable.crop(tmin=5.0) + # the failure left the object exactly as it was + for got, want in zip(variable.get_data(), before_data): + assert_array_equal(got, want) + assert_array_equal(np.array(variable.tmin), before_tmin) + assert_array_equal(np.array(variable.tmax), before_tmax) + assert variable.variable_duration + + +def test_crop_bounds_come_from_retained_samples(variable): + """Test that the stored bounds are sample positions, not the request.""" + # only tmin, so the differing ends keep the object ragged + cropped = variable.copy().crop(tmin=0.013) + assert cropped.variable_duration + # the request fell between samples and was snapped to one + assert not np.isclose(np.atleast_1d(cropped.tmin)[0], 0.013) + for ii in range(len(cropped)): + times = cropped.get_times(ii) + assert times[0] == pytest.approx(np.atleast_1d(cropped.tmin)[ii]) + assert times[-1] == pytest.approx(np.atleast_1d(cropped.tmax)[ii]) + # and the axis still describes the block exactly + assert len(times) == cropped.get_data()[ii].shape[-1] + assert_allclose( + cropped.durations, np.atleast_1d(cropped.tmax) - np.atleast_1d(cropped.tmin) + ) + + +def test_crop_that_equalises_axes_returns_fixed_epochs(variable): + """Test that removing the variation gives an ordinary Epochs back.""" + cropped = variable.copy().crop(tmax=0.5) + assert not cropped.variable_duration + assert isinstance(cropped._data, np.ndarray) + assert cropped._tmin_per_epoch is None + assert cropped._tmax_per_epoch is None + assert isinstance(cropped.tmin, float) + assert isinstance(cropped.tmax, float) + # times is answerable again, and agrees with the data + assert len(cropped.times) == cropped.get_data().shape[-1] + # and the reductions come back on their own, without touching the tables + evoked = cropped.average() + assert evoked.nave == len(cropped) + assert_allclose(evoked.times, cropped.times) + + +def test_crop_keeps_epoch_bookkeeping(variable): + """Test that events, metadata and drop_log travel unchanged.""" + import pandas as pd + + pytest.importorskip("pandas") + variable.metadata = pd.DataFrame(dict(kind=list("abcd"))) + events = variable.events.copy() + drop_log = variable.drop_log + selection = variable.selection.copy() + + cropped = variable.copy().crop(tmin=0.0, tmax=0.4) + assert_array_equal(cropped.events, events) + assert cropped.drop_log == drop_log + assert_array_equal(cropped.selection, selection) + assert list(cropped.metadata["kind"]) == list("abcd") + + +def test_crop_does_not_fall_back_to_as_fixed(variable): + """Test that cropping is native, never a padded copy.""" + + def _boom(*args, **kwargs): + raise AssertionError("crop() fell back to as_fixed()") + + variable.as_fixed = _boom + cropped = variable.crop(tmin=0.0) + assert not np.isnan(np.concatenate(cropped.get_data(), axis=-1)).any() + + +def test_crop_refuses_when_rejection_windows_are_set(variable): + """Test that a stray rejection window is refused, not compared to an array.""" + variable.reject_tmin = 0.0 # the constructor forbids this; be defensive + with pytest.raises(NotImplementedError, match="reject_tmin is not implemented"): + variable.crop(tmin=0.0) From 90aee6df10c8795d386bc7ff3a067615822c77b5 Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 20:40:22 +0200 Subject: [PATCH 06/11] FIX: browser defects found by a corner-case sweep Route variable-duration browsing on a backend capability flag rather than the backend name, so mne-qt-browser can opt in. Duck-typed the way BrowserBase._has_time_slice already is, so an older mne-qt-browser still declines. Six tests whose skip reason claimed matplotlib-only now run on both backends. Three of these fixes are regressions on the *equal-duration* path, not ragged-only: - _recompute_epochs_vlines computed an unclamped sample offset, so a click in an epoch's last half sample rounded one sample past its end, a latency no epoch holds. Every line was dropped while the readout still showed the out-of-range value. The Qt backend kept its old path behind a guard; the matplotlib rewrite had none, so fixed-duration data went through the new code. - _draw_traces rebuilt the visible-epoch list by searchsorting the time range, which drops the last epoch when it holds one sample, and raises when that is the only visible epoch. Ask the view instead. - The colour band mask excluded each epoch's own first sample, so a one-sample epoch was drawn in its neighbour's colour and the window's first sample was never painted at all. Ragged-only: - _create_epoch_histogram called np.ptp(..., axis=2) on a list of arrays. Peak-to-peak is per trial, so compute it per epoch. - _getitem moved the per-epoch bounds but never re-derived the union time axis, so as_fixed() kept padding out to epochs that had been dropped: epochs[0] of a 100-sample epoch returned (1, 3, 280) with 540 NaN and n_contributing == 0 at 180 time points. crop() already re-derives it. - drop_bad(reject=...) reached Epochs.times and raised an internal error with no classification; it now declines clearly. The no-arg call still short-circuits, which _concatenate_epochs relies on. Also silence two ty diagnostics in _crop_variable that predate this branch. Verified against the pre-PR commit across 7,697 recorded states and 377 figures per environment: all 32 fields that determine what the reader sees are bit-identical, and the only movement is the vline landing on a real sample instead of between two. --- mne/epochs.py | 17 +++++++++++++++-- mne/utils/mixin.py | 9 +++++++++ mne/viz/_figure.py | 27 ++++++++++++++++++++++++++- mne/viz/_mpl_figure.py | 22 ++++++++++++++++++---- mne/viz/epochs.py | 12 +++--------- mne/viz/tests/test_epochs.py | 22 +++++++--------------- 6 files changed, 78 insertions(+), 31 deletions(-) diff --git a/mne/epochs.py b/mne/epochs.py index 556dc8c8f1d..40d14937ce8 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -1978,6 +1978,19 @@ def drop_bad( flat = self.flat if any(isinstance(rej, str) and rej != "existing" for rej in (reject, flat)): raise ValueError('reject and flat, if strings, must be "existing"') + if self._variable_duration and (reject or flat): + # the no-arg call has already short-circuited above, so reaching + # here means real thresholds; amplitude rejection would run per + # epoch, but it goes through `times` on the way and there is no + # ragged path for it yet + raise NotImplementedError( + "drop_bad() with reject or flat is not implemented for " + "variable-duration epochs. Amplitude rejection is per-trial " + "and could work here, but the preloaded path needs one shared " + "time axis, which these epochs do not have. Pass reject= to " + "Epochs() at construction instead, which does apply it per " + "epoch. See https://github.com/mne-tools/mne-python/issues/14206." + ) self._reject_setup(reject, flat, allow_callable=True) self._get_data(out=False, verbose=verbose) return self @@ -2855,8 +2868,8 @@ def _crop_variable(self, tmin, tmax, include_tmax): self._tmax_per_epoch = None # ty: ignore[invalid-assignment] start_idx, stop_idx = first_idx[0], first_idx[0] + lengths[0] - 1 else: - start_idx = int(round(starts.min() * sfreq)) - stop_idx = int(round(stops.max() * sfreq)) + start_idx = int(round(float(np.min(starts)) * sfreq)) + stop_idx = int(round(float(np.max(stops)) * sfreq)) self._raw_times = np.arange(start_idx, stop_idx + 1) / sfreq self._set_times(self._raw_times) diff --git a/mne/utils/mixin.py b/mne/utils/mixin.py index 648aefa73c6..46d3f2e205a 100644 --- a/mne/utils/mixin.py +++ b/mne/utils/mixin.py @@ -277,6 +277,15 @@ def _getitem( # an ndarray takes a slice or an index array equally well inst._tmin_per_epoch = inst._tmin_per_epoch[select] inst._tmax_per_epoch = inst._tmax_per_epoch[select] + # the union window is defined by those bounds, so it has to be + # re-derived here as `crop` does; otherwise `as_fixed` keeps + # padding out to epochs that are no longer present + if len(inst._tmin_per_epoch): + sfreq = float(inst.info["sfreq"]) + start_idx = int(round(inst._tmin_per_epoch.min() * sfreq)) + stop_idx = int(round(inst._tmax_per_epoch.max() * sfreq)) + inst._raw_times = np.arange(start_idx, stop_idx + 1) / sfreq + inst._set_times(inst._raw_times) if drop_event_id: # update event id to reflect new content of inst inst.event_id = { diff --git a/mne/viz/_figure.py b/mne/viz/_figure.py index 1e6f324abd4..d1a0899040d 100644 --- a/mne/viz/_figure.py +++ b/mne/viz/_figure.py @@ -710,7 +710,11 @@ def _create_epoch_histogram(self): """Create peak-to-peak histogram of channel amplitudes.""" epochs = self.mne.inst data = OrderedDict() - ptp = np.ptp(epochs.get_data(copy=False), axis=2) + # per epoch, so that variable-duration epochs (a list of arrays) + # work too; peak-to-peak is a per-trial reduction either way + ptp = np.array( + [np.ptp(epoch, axis=-1) for epoch in epochs.get_data(copy=False)] + ) for ch_type in ("eeg", "mag", "grad"): if ch_type in epochs: data[ch_type] = ptp.T[self.mne.ch_types == ch_type].ravel() @@ -808,6 +812,27 @@ def _load_backend(backend_name): return backend +def _check_variable_duration_backend(): + """Raise unless the active browser backend can draw ragged epochs. + + Matplotlib always can. The Qt backend gained the ability in + mne-qt-browser 0.8, which announces it with a module-level flag, so an + older one declines here rather than drawing a wrong picture from a + boundary model it does not know about. + """ + backend_name = get_browser_backend() + if backend_name == "matplotlib": + return + module = _load_backend(backend_name) + if not getattr(module, "_SUPPORTS_VARIABLE_DURATION", False): + raise NotImplementedError( + f"Browsing variable-duration epochs is not implemented for the " + f"{backend_name} backend of this version, only for matplotlib. " + "Upgrade mne-qt-browser, or select matplotlib with " + 'mne.viz.set_browser_backend("matplotlib").' + ) + + def _get_browser(show, block, **kwargs): """Instantiate a new MNE browse-style figure.""" from .utils import _get_figsize_from_config diff --git a/mne/viz/_mpl_figure.py b/mne/viz/_mpl_figure.py index 2ab30e079c9..18151451a6c 100644 --- a/mne/viz/_mpl_figure.py +++ b/mne/viz/_mpl_figure.py @@ -2263,8 +2263,10 @@ def _draw_traces(self): # check for bad epochs time_range = (self.mne.times + self.mne.first_time)[[0, -1]] if self.mne.instance_type == "epochs": - epoch_ix = np.searchsorted(self.mne.boundary_times, time_range) - epoch_ix = np.arange(epoch_ix[0], epoch_ix[1]) + # ask the view directly: deriving this from the time range drops + # the last epoch whenever it holds a single sample, because that + # sample's time is its own left boundary + epoch_ix = np.arange(*self._get_epoch_ix_range()) epoch_nums = self.mne.inst.selection[epoch_ix[0] : epoch_ix[-1] + 1] (visible_bad_epoch_ix,) = np.isin(epoch_nums, self.mne.bad_epochs).nonzero() while len(self.mne.epoch_traces): @@ -2334,7 +2336,11 @@ def _draw_traces(self): _starts = self.mne.boundary_times[epoch_ix][bool_ixs] _stops = self.mne.boundary_times[epoch_ix + 1][bool_ixs] for _start, _stop in zip(_starts, _stops): - _mask = np.logical_and(_start < this_times, this_times <= _stop) + # inclusive at both ends: an epoch owns its own first + # sample, and a one-sample epoch has nothing else + _mask = np.logical_and( + _start <= this_times, this_times <= _stop + ) mask = mask | _mask _times = np.ma.masked_array(this_times, mask=~mask) # always use the existing traces first @@ -2442,7 +2448,15 @@ def _recompute_epochs_vlines(self, xdata): len(boundary_times) - 2, ) ) - offset = round((xdata - boundary_times[clicked_ix]) * sfreq) + # clamp into the clicked epoch: a click in its last half sample would + # otherwise round up to one sample past its end, which no epoch holds, + # and every line would be dropped + n_samp = int( + round((boundary_times[clicked_ix + 1] - boundary_times[clicked_ix]) * sfreq) + ) + offset = int( + np.clip(round((xdata - boundary_times[clicked_ix]) * sfreq), 0, n_samp - 1) + ) latency = self.mne.epoch_tmins[clicked_ix] + offset / sfreq ix_start, ix_stop = self._get_epoch_ix_range() xs = list() diff --git a/mne/viz/epochs.py b/mne/viz/epochs.py index 3f6a4c0f90d..9f182b0f64e 100644 --- a/mne/viz/epochs.py +++ b/mne/viz/epochs.py @@ -1133,15 +1133,9 @@ def plot_epochs( ) if epochs._variable_duration: - from ._figure import get_browser_backend - - backend_name = get_browser_backend() - if backend_name != "matplotlib": - raise NotImplementedError( - f"Browsing variable-duration epochs is not implemented for the " - f"{backend_name} backend yet, only for matplotlib. Select it " - 'with mne.viz.set_browser_backend("matplotlib").' - ) + from ._figure import _check_variable_duration_backend + + _check_variable_duration_backend() fig = _get_browser(show=show, block=block, **params) diff --git a/mne/viz/tests/test_epochs.py b/mne/viz/tests/test_epochs.py index 78fbcf04513..c0d9734fc35 100644 --- a/mne/viz/tests/test_epochs.py +++ b/mne/viz/tests/test_epochs.py @@ -554,8 +554,6 @@ def variable_epochs(): def test_plot_variable_duration_is_native(variable_epochs, browser_backend): """Test that browsing ragged epochs neither warns nor pads.""" - if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") def _boom(*args, **kwargs): raise AssertionError("plot() fell back to as_fixed() instead of browsing") @@ -569,8 +567,6 @@ def _boom(*args, **kwargs): def test_plot_variable_duration_boundaries(variable_epochs, browser_backend): """Test that the browser lays epochs end to end at their true lengths.""" - if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") fig = variable_epochs.plot(n_epochs=2) assert_allclose(fig.mne.boundary_times, _boundaries()) assert_array_equal(fig.mne.boundary_samples, np.r_[0, np.cumsum(LENGTHS_VAR)]) @@ -583,8 +579,6 @@ def test_plot_variable_duration_window_spans_whole_epochs( variable_epochs, browser_backend ): """Test that n_epochs means epochs, not a representative duration.""" - if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") boundaries = _boundaries() for n_epochs in (1, 2, 3, 4): fig = variable_epochs.plot(n_epochs=n_epochs) @@ -594,8 +588,6 @@ def test_plot_variable_duration_window_spans_whole_epochs( def test_plot_variable_duration_data_is_unpadded(variable_epochs, browser_backend): """Test that a view holds exactly the source samples, in order.""" - if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") fig = variable_epochs.plot(n_epochs=2) source = variable_epochs.get_data() for keys in ([], ["right"], ["right", "right"]): @@ -616,8 +608,6 @@ def test_plot_variable_duration_data_is_unpadded(variable_epochs, browser_backen def test_plot_variable_duration_navigation(variable_epochs, browser_backend): """Test that arrow keys move by epochs and land on real boundaries.""" - if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") boundaries = _boundaries() fig = variable_epochs.plot(n_epochs=2) assert fig.mne.t_start == pytest.approx(boundaries[0]) @@ -642,8 +632,6 @@ def test_plot_variable_duration_navigation(variable_epochs, browser_backend): def test_plot_variable_duration_home_end(variable_epochs, browser_backend): """Test that home/end change the epoch count and recompute the duration.""" - if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") boundaries = _boundaries() fig = variable_epochs.plot(n_epochs=2) assert fig.mne.duration == pytest.approx(boundaries[2]) @@ -672,7 +660,7 @@ def test_plot_variable_duration_hscroll_patches(variable_epochs, browser_backend def test_plot_variable_duration_bad_epoch(variable_epochs, browser_backend): """Test that a click finds the right epoch when the widths differ.""" if browser_backend.name != "matplotlib": - pytest.skip("variable-duration browsing is matplotlib-only") + pytest.skip("epoch marking by click is matplotlib-specific") boundaries = _boundaries() fig = variable_epochs.plot(n_epochs=4) y = fig.mne.traces[0].get_ydata()[0] @@ -723,10 +711,14 @@ def test_plot_variable_duration_events(browser_backend): assert_allclose(got, want) -def test_plot_variable_duration_refuses_other_backends(variable_epochs, monkeypatch): - """Test that non-matplotlib backends decline rather than fail obscurely.""" +def test_plot_variable_duration_refuses_old_backends(variable_epochs, monkeypatch): + """Test that a backend without the boundary model declines, not fails.""" import mne.viz._figure + class _OldBackend: # an mne-qt-browser that predates the boundary model + pass + monkeypatch.setattr(mne.viz._figure, "get_browser_backend", lambda: "qt") + monkeypatch.setattr(mne.viz._figure, "_load_backend", lambda name: _OldBackend()) with pytest.raises(NotImplementedError, match="not implemented for the qt"): variable_epochs.plot() From 60656fab43a21852b171175c0f14ff35972ff6f3 Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 22:17:12 +0200 Subject: [PATCH 07/11] MAINT: shorten new private helper docstrings Reduce the private helpers this branch adds to one-line summaries, dropping the Parameters/Returns/Notes sections, and cut the comments that only restated the line below them. Kept as they were: `_get_epoch_ix_range` and `_check_variable_duration_backend` in viz/_figure.py, which record the mne-qt-browser flag contract and the single-source-of-truth invariant behind the shape assertion in `_update_data`. `_decim_slice` in `__init__` was commented as avoiding a densifying `decimate()`; `decimate` is in `_VARIABLE_NOT_IMPLEMENTED` and raises, so the comment now says what the assignment is actually for. No behaviour change: stripping docstrings leaves every touched file with an identical AST, and tests, tutorial and doc/ are untouched. --- mne/channels/channels.py | 6 +- mne/epochs.py | 210 +++++---------------------------------- mne/utils/mixin.py | 4 +- mne/viz/_figure.py | 27 +---- mne/viz/_mpl_figure.py | 8 +- 5 files changed, 31 insertions(+), 224 deletions(-) diff --git a/mne/channels/channels.py b/mne/channels/channels.py index b4e40ad2b71..7ab98ddd261 100644 --- a/mne/channels/channels.py +++ b/mne/channels/channels.py @@ -637,10 +637,8 @@ def _pick_drop_channels(self, idx, *, verbose=None): axis = -2 if hasattr(self, "_data"): # skip non-preloaded Raw if isinstance(self._data, list): - # variable-duration epochs: one array per epoch, channels are - # regular within each, so the pick applies the same way to all. - # Replacing the contents rather than the attribute keeps `_data` - # an ndarray everywhere else this mixin is used. + # replacing the contents rather than the attribute keeps `_data` + # an ndarray everywhere else this mixin is used self._data[:] = [epoch.take(idx, axis=axis) for epoch in self._data] else: self._data = self._data.take(idx, axis=axis) diff --git a/mne/epochs.py b/mne/epochs.py index 40d14937ce8..7968ab654c9 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -387,26 +387,7 @@ def _handle_event_repeated(events, event_id, event_repeated, selection, drop_log def _check_variable_bounds(tmin, tmax, n_events): - """Normalize ``tmin``/``tmax``, which may be given per event. - - Parameters - ---------- - tmin : float | array of float - Start time(s) in seconds. - tmax : float | array of float - End time(s) in seconds. - n_events : int - Number of events, used to check the length of array inputs. - - Returns - ------- - tmin : float | array of float - The validated start time(s). - tmax : float | array of float - The validated end time(s). - variable_duration : bool - Whether either bound was given per event. - """ + """Normalize ``tmin``/``tmax``, which may be given per event.""" arrays = [np.ndim(bound) > 0 for bound in (tmin, tmax)] if not any(arrays): return tmin, tmax, False @@ -442,21 +423,7 @@ def _check_variable_bounds(tmin, tmax, n_events): def _check_variable_data(data, tmin, tmax, events, sfreq): - """Validate per-epoch data for variable-duration epochs. - - Parameters - ---------- - data : list of array - Per-epoch data, each of shape ``(n_channels, n_times_i)``. - tmin : array of float - Per-epoch start times in seconds. - tmax : array of float - Per-epoch end times in seconds. - events : array of int - The events. - sfreq : float - The sampling frequency in Hz. - """ + """Validate per-epoch data for variable-duration epochs.""" if isinstance(data, np.ndarray) and data.ndim == 3: data = list(data) if not isinstance(data, list | tuple): @@ -491,21 +458,7 @@ def _check_variable_data(data, tmin, tmax, events, sfreq): def _check_variable_unsupported(*, baseline, reject_tmin, reject_tmax, decim, preload): - """Reject options this implementation does not yet handle. - - Parameters - ---------- - baseline : tuple | None - The requested baseline. - reject_tmin : float | None - Start of the rejection window. - reject_tmax : float | None - End of the rejection window. - decim : int - Decimation factor. - preload : bool - Whether the data is preloaded. - """ + """Reject options this implementation does not yet handle.""" if baseline is not None: raise NotImplementedError( "Baseline correction is not implemented for variable-duration " @@ -531,18 +484,7 @@ def _check_variable_unsupported(*, baseline, reject_tmin, reject_tmax, decim, pr def _is_variable_duration_data(data): - """Whether ``data`` is a sequence of arrays with differing lengths. - - Parameters - ---------- - data : array | list of array - Candidate epoch data. - - Returns - ------- - variable : bool - True if the entries cannot share one time axis. - """ + """Whether ``data`` is a sequence of arrays with differing lengths.""" if isinstance(data, np.ndarray): return False if not isinstance(data, list | tuple) or len(data) == 0: @@ -777,8 +719,7 @@ def __init__( self.metadata = metadata # do not set self.events here, let subclass do it - # Variable-duration epochs: tmin and/or tmax may be given per event. The - # scalar path below is untouched; only _variable_duration is new. + # Variable-duration epochs: tmin and/or tmax may be given per event. tmin, tmax, self._variable_duration = _check_variable_bounds( tmin, tmax, @@ -891,7 +832,8 @@ def __init__( if not self._variable_duration: self.decimate(decim) else: - # decim != 1 is rejected above; decimate() would densify _data + # decim != 1 is rejected above, and decimate() is wrapped to raise + # for these, so set the slice it would otherwise have set self._decim_slice = slice(None, None, None) # baseline correction: replace `None` tuple elements with actual times @@ -1132,24 +1074,7 @@ def as_fixed(self, pad_value=np.nan): return out, n_contributing def _get_variable_data(self, *, picks=None, item=None, copy=True): - """Return per-epoch data when durations vary. - - Parameters - ---------- - picks : str | array-like | slice | None - Channels to include. - item : slice | array-like | str | list | None - Epochs to include. - copy : bool - Whether to copy the data. - - Returns - ------- - data : list of array - One ``(n_channels, n_times_i)`` array per epoch. A single array is - not returned, because there is no length every epoch shares; use - :meth:`as_fixed` to obtain one. - """ + """Return per-epoch data when durations vary, as one array per epoch.""" if item is None: item = slice(None) sel = np.arange(len(self.events))[item] if not isinstance(item, str) else None @@ -1181,7 +1106,6 @@ def load_data(self) -> Self: if self.preload: return self if self._variable_duration: - # a list, one array per epoch; see the note on _data above self._data = self._load_variable_from_raw() else: self._data = self._get_data() @@ -1352,18 +1276,7 @@ def _reject_setup(self, reject, flat, *, allow_callable=False): self._reject_time = slice(reject_imin, reject_imax) def _load_variable_from_raw(self): - """Read every epoch at its own length, dropping bad ones. - - Returns - ------- - data : list of array - One ``(n_channels, n_times_i)`` array per retained epoch. - - Notes - ----- - Mirrors the drop bookkeeping of :meth:`drop_bad` for the fixed-duration - path, but collects a list because there is no length the epochs share. - """ + """Read every epoch at its own length, mirroring ``drop_bad`` bookkeeping.""" detrend_picks = self._detrend_picks drop_log = list(self.drop_log) good_idx, out = [], [] @@ -1397,18 +1310,7 @@ def _load_variable_from_raw(self): return out def _n_times_per_epoch(self, idx): - """Return the number of samples in one epoch. - - Parameters - ---------- - idx : int - Index of the epoch. - - Returns - ------- - n_times : int - Number of samples. - """ + """Return the number of samples in one epoch.""" if not self._variable_duration: return len(self.times) sfreq = float(self.info["sfreq"]) @@ -1980,9 +1882,7 @@ def drop_bad( raise ValueError('reject and flat, if strings, must be "existing"') if self._variable_duration and (reject or flat): # the no-arg call has already short-circuited above, so reaching - # here means real thresholds; amplitude rejection would run per - # epoch, but it goes through `times` on the way and there is no - # ragged path for it yet + # here means real thresholds raise NotImplementedError( "drop_bad() with reject or flat is not implemented for " "variable-duration epochs. Amplitude rejection is per-trial " @@ -2765,23 +2665,7 @@ def crop( return self def _crop_variable(self, tmin, tmax, include_tmax): - """Crop each epoch on its own time axis. - - The requested window is a physical interval in seconds, so it is applied - to every epoch independently and clamped to the epoch's own bounds when - it reaches past them. Nothing is padded, interpolated or aligned; an - epoch that the window misses entirely makes the whole call fail, the - same way it would for that epoch on its own. - - Parameters - ---------- - tmin : float | None - Start of the window, or ``None`` for each epoch's own start. - tmax : float | None - End of the window, or ``None`` for each epoch's own end. - include_tmax : bool - Whether to keep the sample at ``tmax``. - """ + """Crop each epoch on its own time axis.""" for name in ("reject_tmin", "reject_tmax"): if getattr(self, name, None) is not None: raise NotImplementedError( @@ -2839,7 +2723,6 @@ def _crop_variable(self, tmin, tmax, include_tmax): "each of those epochs' own last sample." ) - # Second pass: apply ragged = self._data assert ragged is not None # variable-duration epochs are always preloaded starts = np.empty(len(masks)) @@ -4169,16 +4052,16 @@ def _events_from_annotations(raw, events, event_id, annotations, on_missing): return events, event_id, annotations -#: Methods whose result is only looked at. They warn and run on ``as_fixed()``, -#: which is enough for inspection because the padding is visible to whoever is -#: looking. Anything numeric is not in this table: see the two below. +# Methods whose result is only looked at. They warn and run on ``as_fixed()``, +# which is enough for inspection because the padding is visible to whoever is +# looking. Anything numeric is not in this table: see the two below. _VARIABLE_FALLBACK = { "to_data_frame": "", } -#: Methods that combine epochs across a shared time axis. Padding them makes the -#: number of contributing epochs a function of time, which no scalar ``nave`` can -#: describe, so they ask for a policy instead of inventing one. See :gh:`14206`. +# Methods that combine epochs across a shared time axis. Padding them makes the +# number of contributing epochs a function of time, which no scalar ``nave`` can +# describe, so they ask for a policy instead of inventing one. See :gh:`14206`. _VARIABLE_NEEDS_POLICY = { "average": "averaging", "standard_error": "estimating the standard error", @@ -4188,9 +4071,9 @@ def _events_from_annotations(raw, events, event_id, annotations, on_missing): "compute_psd": "computing a spectrum", } -#: Methods that are mathematically per-trial and simply have no ragged -#: implementation yet. Running them on a padded copy would return a wrong answer -#: rather than a slow one, so they raise until implemented natively. +# Methods that are mathematically per-trial and simply have no ragged +# implementation yet. Running them on a padded copy would return a wrong answer +# rather than a slow one, so they raise until implemented natively. _VARIABLE_NOT_IMPLEMENTED = { "filter": "filtering", "apply_function": "applying a function", @@ -4207,22 +4090,7 @@ def _events_from_annotations(raw, events, event_id, annotations, on_missing): def _wrap_variable_fallback(func, name, note): - """Warn and fall back to ``as_fixed()`` for variable-duration epochs. - - Parameters - ---------- - func : callable - The original method. - name : str - Its name, used to look it up on the fixed-duration copy. - note : str - Extra sentence appended to the warning, or an empty string. - - Returns - ------- - wrapper : callable - The wrapped method. - """ + """Warn and fall back to ``as_fixed()`` for variable-duration epochs.""" @wraps(func) def wrapper(self, *args, **kwargs): @@ -4244,22 +4112,7 @@ def wrapper(self, *args, **kwargs): def _raise_needs_policy(func, name, what): - """Raise for reductions across a time axis the epochs do not share. - - Parameters - ---------- - func : callable - The original method. - name : str - Its name. - what : str - Short description of the operation, used in the message. - - Returns - ------- - wrapper : callable - The wrapped method. - """ + """Raise for reductions across a time axis the epochs do not share.""" @wraps(func) def wrapper(self, *args, **kwargs): @@ -4279,22 +4132,7 @@ def wrapper(self, *args, **kwargs): def _raise_not_implemented(func, name, what): - """Raise for per-trial operations with no ragged implementation yet. - - Parameters - ---------- - func : callable - The original method. - name : str - Its name. - what : str - Short description of the operation, used in the message. - - Returns - ------- - wrapper : callable - The wrapped method. - """ + """Raise for per-trial operations with no ragged implementation yet.""" @wraps(func) def wrapper(self, *args, **kwargs): diff --git a/mne/utils/mixin.py b/mne/utils/mixin.py index 46d3f2e205a..b276f524c5a 100644 --- a/mne/utils/mixin.py +++ b/mne/utils/mixin.py @@ -223,8 +223,7 @@ def _getitem( inst = self.copy() if copy else self if self._data is not None: if isinstance(self._data, list): - # variable-duration epochs hold one array per epoch, so there is - # nothing to copy into; copy() already produced the list + # one array per epoch, so there is nothing to copy into inst._data = [d.copy() for d in self._data] if copy else self._data else: np.copyto(inst._data, self._data, casting="no") @@ -272,7 +271,6 @@ def _getitem( # ensure that each Epochs instance owns its own data so we can # resize later if necessary inst._data = np.require(inst._data[select], requirements=["O"]) - # per-event bounds travel with the epochs they describe if getattr(inst, "_variable_duration", False): # an ndarray takes a slice or an index array equally well inst._tmin_per_epoch = inst._tmin_per_epoch[select] diff --git a/mne/viz/_figure.py b/mne/viz/_figure.py index d1a0899040d..2cbe146f715 100644 --- a/mne/viz/_figure.py +++ b/mne/viz/_figure.py @@ -45,29 +45,7 @@ def __init__(self, **kwargs): def _epoch_window(boundary_times, start_ix, n_epochs): - """Return the start time and duration of a window of whole epochs. - - Epochs may differ in duration, so a window of ``n_epochs`` of them spans - whatever lies between the two boundaries rather than a fixed number of - seconds. ``start_ix`` is clamped so the requested epochs stay visible when - the object is long enough to allow it. - - Parameters - ---------- - boundary_times : array - Cumulative epoch edges in seconds, including both ends. - start_ix : int - Index of the first epoch to show. - n_epochs : int - Number of epochs to show. - - Returns - ------- - t_start : float - Time of the first boundary. - duration : float - Seconds spanned by the requested epochs. - """ + """Return the start time and duration of a window of whole epochs.""" n_total = len(boundary_times) - 1 n_epochs = int(np.clip(n_epochs, 1, n_total)) start_ix = int(np.clip(start_ix, 0, n_total - n_epochs)) @@ -412,8 +390,7 @@ def _get_start_stop(self): # update time start_sec = self.mne.t_start - self.mne.first_time if self.mne.is_epochs: - # take the samples the visible epochs really hold, so that this - # agrees with _load_data by construction rather than by arithmetic + # this agrees with _load_data by construction, not by arithmetic ix_start, ix_stop = self._get_epoch_ix_range() start = int(self.mne.boundary_samples[ix_start]) stop = int(self.mne.boundary_samples[ix_stop]) diff --git a/mne/viz/_mpl_figure.py b/mne/viz/_mpl_figure.py index 18151451a6c..3ad5989c80d 100644 --- a/mne/viz/_mpl_figure.py +++ b/mne/viz/_mpl_figure.py @@ -791,8 +791,7 @@ def _keypress(self, event): old_t_start = self.mne.t_start direction = 1 if key.endswith("right") else -1 if self.mne.is_epochs: - # step whole epochs, since they need not share a duration: one - # epoch normally, a whole window with shift + # step whole epochs, since they need not share a duration step = self.mne.n_epochs if key.startswith("shift") else 1 ix_start, _ = self._get_epoch_ix_range() self.mne.t_start, self.mne.duration = _epoch_window( @@ -841,8 +840,7 @@ def _keypress(self, event): self.mne.n_epochs = int( np.clip(self.mne.n_epochs + dur_delta, 1, len(self.mne.inst)) ) - # the epochs added or removed have their own durations, so ask - # the boundaries how many seconds that actually is + # the epochs added or removed have their own durations self.mne.t_start, new_dur = _epoch_window( self.mne.boundary_times, ix_start, self.mne.n_epochs ) @@ -1012,7 +1010,6 @@ def _mouse_move(self, event): time = np.clip(time, self.mne.first_time, max_time) if self.mne.is_epochs: ix = np.searchsorted(self.mne.boundary_times[1:], time, side="right") - # the epochs from here on have their own durations time, self.mne.duration = _epoch_window( self.mne.boundary_times, ix, self.mne.n_epochs ) @@ -1934,7 +1931,6 @@ def _check_update_hscroll_clicked(self, event): time = np.clip(time, self.mne.first_time, max_time) if self.mne.is_epochs: ix = np.searchsorted(self.mne.boundary_times[1:], time, side="right") - # the epochs from here on have their own durations time, self.mne.duration = _epoch_window( self.mne.boundary_times, ix, self.mne.n_epochs ) From e3934a20320725d652f9c3165dff789c69a2641d Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 22:22:05 +0200 Subject: [PATCH 08/11] FIX: refuse get_data arguments that ragged epochs cannot honour `get_data()` dispatched to `_get_variable_data(picks, item, copy)` for variable-duration epochs, so `units`, `tmin` and `tmax` were accepted and then discarded: `get_data(units="uV")` returned volts and `get_data(tmin=..., tmax=...)` returned whole epochs. Raise instead, naming the argument and pointing at `as_fixed()`, which supports all three. The docstring promised a 3D array on both paths; it and the return annotation now cover the list of one array per epoch that ragged epochs return. `save()` asserts the array case it already guarantees, since it refuses ragged epochs. --- mne/epochs.py | 18 +++++++++++++++--- mne/tests/test_epochs_variable_duration.py | 17 +++++++++++++++++ 2 files changed, 32 insertions(+), 3 deletions(-) diff --git a/mne/epochs.py b/mne/epochs.py index 7968ab654c9..a322b98cfc4 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -2380,9 +2380,12 @@ def get_data( *, copy: bool = True, verbose: bool | str | int | None = None, - ) -> np.ndarray: + ) -> np.ndarray | list[np.ndarray]: """Get all epochs as a 3D array. + Epochs of different durations have no shared time axis, so those are + returned as one array per epoch instead; see :meth:`as_fixed`. + Parameters ---------- %(picks_all)s @@ -2425,11 +2428,19 @@ def get_data( Returns ------- - data : array of shape (n_epochs, n_channels, n_times) + data : array of shape (n_epochs, n_channels, n_times) | list of array The epochs data. Will be a copy when ``copy=True`` and will be a view - when possible when ``copy=False``. + when possible when ``copy=False``. When durations vary, one + ``(n_channels, n_times_i)`` array per epoch. """ if self._variable_duration: + for name, value in (("units", units), ("tmin", tmin), ("tmax", tmax)): + if value is not None: + raise NotImplementedError( + f"get_data() with {name} is not implemented for " + "variable-duration epochs; it would have to be applied " + "per epoch. Call as_fixed() first to get one array." + ) return self._get_variable_data(picks=picks, item=item, copy=copy) return self._get_data( picks=picks, item=item, units=units, tmin=tmin, tmax=tmax, copy=copy @@ -2863,6 +2874,7 @@ def save( total_size = 0 else: d = self[0].get_data(copy=False) + assert isinstance(d, np.ndarray) # save() refuses ragged epochs # this should be guaranteed by subclasses assert d.dtype in (">f8", "c16", " Date: Thu, 27 Aug 2026 22:27:43 +0200 Subject: [PATCH 09/11] MAINT: two follow-ups to the docstring pass `_load_variable_from_raw` returns a list, which its one-line summary now says, since the call site comment that said so went with the sections. `# First pass:` labelled a pair whose second half was a bare `# Second pass: apply` and was removed; drop the label rather than restore it. --- mne/epochs.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/mne/epochs.py b/mne/epochs.py index a322b98cfc4..e9b5f0a6931 100644 --- a/mne/epochs.py +++ b/mne/epochs.py @@ -1276,7 +1276,7 @@ def _reject_setup(self, reject, flat, *, allow_callable=False): self._reject_time = slice(reject_imin, reject_imax) def _load_variable_from_raw(self): - """Read every epoch at its own length, mirroring ``drop_bad`` bookkeeping.""" + """Read every epoch at its own length into a list, mirroring ``drop_bad``.""" detrend_picks = self._detrend_picks drop_log = list(self.drop_log) good_idx, out = [], [] @@ -2686,8 +2686,8 @@ def _crop_variable(self, tmin, tmax, include_tmax): ) sfreq = float(self.info["sfreq"]) - # First pass: work out every selection while changing nothing, so that a - # window that misses one epoch leaves the object as it was. + # Work out every selection before changing anything, so a window that + # misses one epoch leaves the object as it was. masks = list() clamped_tmin = clamped_tmax = False for ii in range(len(self.events)): From 75d6e5aa8b9e3b69523ffbbb3c692a274e56aa48 Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 22:41:24 +0200 Subject: [PATCH 10/11] FIX: skip the metadata crop test when pandas is missing `test_crop_keeps_epoch_bookkeeping` called `pytest.importorskip("pandas")` after `import pandas as pd`, so the hard import raised first and the guard never ran. This failed the minimal build. --- mne/tests/test_epochs_variable_duration.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/mne/tests/test_epochs_variable_duration.py b/mne/tests/test_epochs_variable_duration.py index 667f5f4088c..347b24b9686 100644 --- a/mne/tests/test_epochs_variable_duration.py +++ b/mne/tests/test_epochs_variable_duration.py @@ -670,9 +670,7 @@ def test_crop_that_equalises_axes_returns_fixed_epochs(variable): def test_crop_keeps_epoch_bookkeeping(variable): """Test that events, metadata and drop_log travel unchanged.""" - import pandas as pd - - pytest.importorskip("pandas") + pd = pytest.importorskip("pandas") variable.metadata = pd.DataFrame(dict(kind=list("abcd"))) events = variable.events.copy() drop_log = variable.drop_log From 1c9e29da73d3e6e437a9f159312a6e12a2b8d3ea Mon Sep 17 00:00:00 2001 From: snesmaeili Date: Thu, 27 Aug 2026 23:05:17 +0200 Subject: [PATCH 11/11] FIX: skip the browser tests when the backend declines ragged epochs The variable-duration browser tests ran under both backends, but `plot()` raises for a qt backend that does not announce `_SUPPORTS_VARIABLE_DURATION`, so six of them failed the Ultraslow_PG build. They passed locally only because mne-tools/mne-qt-browser#452 was installed. Ask the guard rather than the backend name, so they run wherever the backend really does support ragged epochs and skip with its own message where it does not. `test_plot_variable_duration_refuses_old_backends` builds its epochs directly, since it must still run when the fixture would skip. --- mne/viz/tests/test_epochs.py | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/mne/viz/tests/test_epochs.py b/mne/viz/tests/test_epochs.py index c0d9734fc35..dd7800d271a 100644 --- a/mne/viz/tests/test_epochs.py +++ b/mne/viz/tests/test_epochs.py @@ -547,8 +547,14 @@ def _boundaries(): @pytest.fixture -def variable_epochs(): +def variable_epochs(browser_backend): """Epochs of unequal duration sharing tmin=0.""" + from mne.viz._figure import _check_variable_duration_backend + + try: # skips qt until an mne-qt-browser that announces support is released + _check_variable_duration_backend() + except NotImplementedError as exc: + pytest.skip(str(exc)) return _variable_epochs() @@ -711,7 +717,7 @@ def test_plot_variable_duration_events(browser_backend): assert_allclose(got, want) -def test_plot_variable_duration_refuses_old_backends(variable_epochs, monkeypatch): +def test_plot_variable_duration_refuses_old_backends(monkeypatch): """Test that a backend without the boundary model declines, not fails.""" import mne.viz._figure @@ -721,4 +727,4 @@ class _OldBackend: # an mne-qt-browser that predates the boundary model monkeypatch.setattr(mne.viz._figure, "get_browser_backend", lambda: "qt") monkeypatch.setattr(mne.viz._figure, "_load_backend", lambda name: _OldBackend()) with pytest.raises(NotImplementedError, match="not implemented for the qt"): - variable_epochs.plot() + _variable_epochs().plot()