ENH: allow Epochs to hold trials of different duration by snesmaeili · Pull Request #14210 · mne-tools/mne-python

snesmaeili · 2026-08-25T07:58:38Z

Draft, following the design discussed in #14206.

tmin and tmax accept (n_events,) arrays. Bounds that carry no actual
variation collapse back to the scalar path, so existing behaviour is unchanged.
EpochsArray also takes a list of (n_channels, n_times_i) arrays and derives
each epoch's tmax from its own length.

New: durations, get_times(epoch), variable_duration, as_fixed().
as_fixed() returns a padded EpochsArray together with the number of epochs
contributing at each time point.

times raises when durations vary, rather than returning the union. Returning it
would leave len(epochs.times) == data.shape[-1] false while looking ordinary.
The union is still available as as_fixed().times.

Method behaviour. Three groups rather than a blanket fallback:

native: pick, drop, __getitem__, shift_time, plot
warn and run on as_fixed(): to_data_frame
raise: the six reductions (average, standard_error, subtract_evoked,
iter_evoked, compute_tfr, compute_psd), plus per-trial operations with no
ragged implementation yet (filter, apply_function, apply_baseline,
crop, decimate, resample, plot_image, plot_topo_image, save,
export)

I originally had the reductions fall back with a warning, as suggested in the
issue. Measuring it changed my mind. On 43 epochs spanning 2.0–3.6 s,
average() returns an Evoked that is 44% NaN from the first drop-out onward,
while nave reports 43 where 3 epochs remain — one short epoch takes out the
whole time point. compute_tfr was worse: it padded and then transformed, which
is the reverse of the order argued for in the issue.

Browsing. Following @drammock's point that Epochs.plot already draws a
pseudo-continuous strip, plot() is native rather than a padded fallback. The
browser's x axis is built from the samples the epochs really hold:

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. A window of k
epochs starting at i spans boundary_times[i + k] - boundary_times[i];
_get_start_stop and _load_data share one index range, so the sample bounds
and the concatenated array agree by construction rather than by arithmetic.
Arrow keys step whole epochs, the scrollbar draws each epoch at its own width,
and a vertical line marks a latency relative to each epoch's own event, omitted
from epochs too short to reach it.

Validated on synthetic epochs of 100, 250, 75 and 180 samples: boundaries are
the cumulative real sample counts, n_times is their sum rather than
n_epochs × max, the loaded window equals np.concatenate of the source blocks
byte for byte, no NaN appears anywhere, and as_fixed() is never called.

Browsing is Matplotlib-only for now — other backends decline with a message
naming it. The PyQtGraph companion is mne-tools/mne-qt-browser#452, which
consumes the same boundary_times / boundary_samples / n_times this puts in
the browser params; relaxing the guard here is a follow-up once that is released.

What was validated. Construction and extraction, plus a per-epoch TFR
pipeline built on get_data(). On all 24 ds004505 subjects (29,546 swing
cycles), every epoch is byte-identical to the raw slice it came from, and
re-running an existing ERSP analysis through the container reproduces the
previously computed maps at 0.000e+00 dB with matching retained counts. Scripts
and per-subject reports:
https://github.com/snesmaeili/meta-mne-python-sprint/tree/main/validation

Tutorial. tutorials/epochs/70_variable_duration_epochs.py builds epochs
from the Sleep Physionet hypnogram durations, so no new dataset is needed and
tools/circleci_download.sh already prefetches it. It shows the container, the
operations that are unaffected, the ones that refuse, browsing bouts of 120, 30,
150, 60 and 90 s at their real widths, and the contributing-count curve from
as_fixed(). Its closing section points at tut-sleep-stage-classif, where
fixed 30 s windows are the right representation, so the two are not read as
competing.

This does not validate any padded path — the methods that would need one raise.

mne/tests/test_epochs.py passes unchanged; 336 tests pass across epochs,
variable-duration and browser tests.

Open questions are in #14206: whether tmin/tmax should stay scalar with
per-epoch bounds under separate names, what the reductions should eventually do
about a contributing count that varies over time, and the FIF representation.

AI assistance: I designed the approach and ran the analyses it is validated
against; Claude Opus 5 wrote the implementation and the tests under my
direction, which I reviewed and tested.

larsoner

Next step I would maybe add a new tutorial so we can see it working. If we need to add a new dataset we could, but maybe better would be to use openneuro-py in the example to download 1 subject's data and process it? MNE-BIDS does something like this and it seems to be okay. (Eventually we'll want to have CircleCI do this ahead of time so triage based on example content and modify _download_all_example_data, but we can do those steps later.)

larsoner · 2026-08-25T08:16:42Z

        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,


Not sure why this would need to be string?

Suggested change

tmin: "float | np.ndarray" = -0.2,

tmin: float | np.ndarray = -0.2,

larsoner · 2026-08-25T08:19:17Z

+
+
 @fill_doc
 class Epochs(BaseEpochs):


I see a lot of decorator mechanics etc. Would things get simpler if we added a EpochsRagged (or some better name) class instead of expanding Epochs itself? I think the decorator idea might have been mine but not sure if it's better or worse than a separate class... I'm thinking this isn't so bad, just want to make sure we thought about another option.

I also spiked EpochsRagged(BaseEpochs) to check whether a separate class would actually simplify this. It doesn't avoid the two shared-module changes: with mixin.py and channels.py restored to upstream, epochs[0:2] fails in _getitem because it assumes ndarray storage, and pick() fails in _pick_drop_channels for the same reason. Avoiding those small shared branches would mean overriding the methods in the subclass and duplicating their selection / drop_log / metadata bookkeeping. Most of the remaining branches in epochs.py are construction and per-epoch-bound handling that would move into the subclass rather than disappear. So at the moment I think keeping one class is simpler, but I'm happy to switch if you prefer the stronger type separation.

The decorator mechanics are independent of that choice. I can replace the dynamic setattr wrappers with explicit variable-duration checks in the affected methods so the behaviour is visible where someone reading average(), filter(), etc. would expect to find it.

For the tutorial, sleep_physionet looks like a good fit and does not require a new dataset. Its hypnogram annotations already carry durations. The existing sleep tutorial intentionally uses chunk_duration=30. for the sleep-stage classification example; on SC4001 that converts 141 annotated sleep-stage bouts spanning 30–1890 s into 653 fixed 30-s events. That gives us a compact way to show the distinction between preserving annotation durations as variable-length epochs and explicitly converting them to a fixed-window representation. I can build the tutorial around one subject and show durations, get_times(), indexing, and as_fixed() with the contributing-count curve

Okay sounds good to me! I think the decorators are simple enough and better than a bunch of repeated checks.

For the new example, CircleCI treats all warnings as errors so you probably need a verbose="error" (easiest way to suppress the warning) during raw read, see:

https://app.circleci.com/pipelines/github/mne-tools/mne-python/33859/workflows/f50bc47a-8410-4a28-8134-67f45e86718e/jobs/86872

(our CIs are under a heavy load so I'm going to kill the currently running ones to save some cycles under the assumption you'll fix this soon-ish and push!)

drammock · 2026-08-25T14:48:44Z

warn and run on as_fixed(): plot, to_data_frame

Just a quick note to say that plot shouldn't be too hard hopefully; our current epochs plotting code effectively concatenates the epochs back into a continuous (raw-like) structure of n_channels x n_times and plots with vertical lines between them. Conceptually nothing changes for ragged epochs (but IDK how hard it will be to make the implementation work, e.g. the scrolling logic may need to be more raw-like than epochs-like...)

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 mne-toolsgh-14206.

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.

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.

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.

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.

larsoner

I see the NOT_IMPLEMENTED list has shrunk, which is good! But it will make it harder to review 😓

I would suggest to stop at the current list, see if there are candidates for simplification, get this plus the mne-qt-browser bit working well, then I find the time to read test manually and merge.

In the meantime, since this implements browsing, if you have more LLM cycles, can you ask a fresh agent to try to find corner cases across multiple interactive use cases (can use qtbot and/or QTest), clicking around, setting channel counts, having lots of epochs, few epochs, vastly different durations, dropped vs not, etc.? This would help ensure that the code here is robust. Might make your laptop kind of unusable for a bit but it will be able to iterate much faster than me...

larsoner · 2026-08-27T14:28:04Z



+def _check_variable_bounds(tmin, tmax, n_events):
+    """Normalize ``tmin``/``tmax``, which may be given per event.


The diff is quite big at +2,030 -109, so I'm hoping / looking for some way to make the diff smaller. I don't see a lot of ways, but one way would actually be to remove the docstrings from these private helpers (other than the first line). LLMs are great at generating a lot of content like this but it ends up needing to be checked and maintained by humans, so if the names and likely types and shapes etc. are already unambiguous from the surrounding context, we have been tending to omit them nowadays. They tend to go out of date quite quickly as well, since private function docstrings are not checked by automated tooling, and it's too easy to forget to update them.

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.

Harness, per-slice reports, screenshots and the draft reply for the sweep behind mne-tools/mne-python#14210 and mne-tools/mne-qt-browser#452. The harness builds every expectation from the source arrays rather than from the object under test, and the fuzzer carries a mutation self-test, so a clean result means something. Raw run logs and caches are ignored; the reports carry the numbers. Findings are triaged against the equal-duration path throughout: anything that also happens without ragged input is recorded as pre-existing rather than fixed. Three defects turned out to regress equal-duration browsing.

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.

`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.

`_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.

`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.

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.

snesmaeili · 2026-08-28T08:26:06Z

Thanks @drammock for your in-person input on this — I agree with the concern about being conservative here.
After digging further into the literature and into the different use cases that have motivated this feature over the years, I think the main thing we need to avoid is treating "variable-duration epochs" as one statistical problem.
The representation problem is relatively clear: sometimes we genuinely need to retain segments with different numbers of samples without truncating them or inventing samples. But what should happen after that — averaging, TFR, PSD, covariance, ICA, source analysis, statistics, etc. — depends heavily on why the durations differ.
So I think the safest direction is:

First support a truthful, explicit representation of variable-domain segments; then make downstream methods ragged-aware only when we can state exactly what they estimate and under what assumptions.
I also agree that an explicitly opt-in / separate variable-duration epoching path may be safer than changing the expectations around ordinary fixed-duration Epochs. The important part for me is that existing fixed-duration workflows remain unchanged.
The main cases we need to distinguish
Case Typical example Why length varies Main analysis family
True incomplete / truncated observation recording dropout, missing tail underlying process continues but part is unobserved missing-data / partially observed functional-data methods
Naturally bounded event spindle, seizure, vocal burst event actually terminates onset/offset locking, duration/event models
Response-terminated interval arithmetic, decision making, self-paced task duration is reaction/processing time stimulus- and response-locked analyses, RT modelling
Repeated process with homologous landmarks gait cycle, respiratory cycle, rhythmic movement same process unfolds at different rates landmark / phase registration
Variable external stimulus sentence, speech segment, movie event stimulus itself has variable duration/content event anchors, encoding/TRF models, variable-domain models
Dense overlapping events fixations, words, tone sequences inter-event intervals vary and responses overlap regression / deconvolution
Variable state sequence sleep architecture, anesthesia states state order and dwell times vary HMM / HSMM / semi-Markov models
Hierarchical multi-stage trial cue → preparation → movement → feedback several segment types, each can have variable duration segment-specific / hierarchical modelling
Those cases can all benefit from the same storage capability, but they should not all get the same implementation of average(), compute_tfr(), etc.
For example:
In gait, a meaningful phase coordinate exists, so landmark registration can be appropriate.
For a sleep spindle, the event ending is not missing data; onset, offset, morphology and duration may themselves be quantities of interest.
For a self-paced arithmetic task, late physical-time samples increasingly come only from slow-response trials.
For whole-night sleep architecture, state bouts are probably better represented and analysed as a state/dwell-time sequence rather than something to average as an Evoked.
For fixations, words or closely spaced tones, the bigger problem may be overlapping neural responses rather than unequal array lengths.
For genuine recording truncation, partially observed functional-data methods may be relevant because an underlying trajectory exists but is unobserved.
So the same representation does not imply the same statistical treatment.
Proposed roadmap
I think we should divide functionality into three groups.

Safe / representation-level functionality
These do not make a cross-trial statistical decision and seem appropriate for an initial implementation.
Area Proposed status
Explicit variable-duration extraction/construction safe
Per-epoch tmin / tmax / duration safe
Per-epoch time vectors such as get_times(epoch) safe
get_data() returning the real per-epoch arrays safe
indexing / slicing / condition selection safe
dropping epochs safe
channel pick/drop/reorder safe
metadata / annotations safe
projections acting only on the spatial/channel axis safe
browsing variable-duration epochs safe
explicit conversion to a rectangular representation safe only when explicitly requested and contributor/mask information is retained
The key constraint should be:

No numerical method should silently convert ragged data into padded data.
If the user explicitly calls something like as_fixed(), that is different: they have requested a representation change and can be given contributor-count / validity information alongside it.
For equal-duration inputs, everything should reduce exactly to current fixed-duration Epochs behaviour.
I think this is approximately where the current PR should stop.

Per-segment methods that are mathematically possible but need careful implementation
These do not necessarily require a shared time axis, so I think they are good candidates for small follow-up PRs.
Method/family Main issue to validate
filter boundary effects depend on segment duration; filtering Raw before segmentation should remain preferred
crop naturally works per segment but changes available support
resample preserves physical duration; should not be confused with time normalization
decimate naturally per segment
Hilbert transform per-segment operation, but short segments and FFT/padding details need parity testing
apply_function naturally maps over individual segments
spatial projections / re-referencing time length itself is irrelevant
baseline correction baseline anchor/window must be defined and available for every segment
artifact rejection rejection probability can depend on epoch duration
source-operator application mathematically sample-wise once inverse/beamformer filters already exist
For these, I think the right acceptance criterion is:

Applying the ragged implementation to epoch i must reproduce applying the current MNE implementation independently to that same epoch.
That gives us a clean numerical parity test.
Artifact rejection deserves special attention
For example, peak-to-peak rejection effectively uses a statistic such as
max(X_i) - min(X_i).
A longer segment has more opportunities to contain an extreme value than a shorter segment. Therefore the probability of rejection can increase with duration even if the underlying noise process is identical.
So rejection is technically per-epoch, but it is not necessarily duration-neutral.
We may eventually want alternatives such as:
fixed-duration quality windows;
proportion of bad samples;
artifact rate per unit time;
or bad-time masks instead of dropping an entire long segment.
That needs more work before we call rejection scientifically equivalent across lengths.

Methods that need a statistical contract before implementation
I would keep these disabled/raising for now.
Method/family Why it is unresolved
average() must define common physical support vs available-case mean vs registered coordinate
standard_error() effective sample count can vary over time
Evoked creation one scalar nave cannot represent time-varying contributor counts
iter_evoked() each input can have a different time axis
subtract_evoked() requires a well-defined common coordinate
averaged TFR requires temporal alignment / reduction semantics
PSD aggregation event-weighted vs time/duration-weighted estimates differ
covariance sample-weighted and epoch-weighted estimates differ
ICA on selected ragged segments longer segments provide more samples unless weighting is deliberately changed
CSP / SPoC covariance estimation and weighting need a policy
time-resolved decoding requires homologous time or phase coordinates
temporal generalization train/test time coordinates must be defined
cluster statistics requires a common test lattice / adjacency
grand averaging requires aligned comparable outputs
dipole fitting from averaged data depends on how the average and covariance were defined
inverse/beamformer estimation covariance and noise-model assumptions need to be settled
connectivity averaged across segments duration, observations and frequency resolution affect estimates
save/export serialization can wait until the in-memory semantics are stable
The important distinction here is:

Some of these are technically easy to make execute, but scientifically easy to make wrong.
Averaging is probably the clearest example
Suppose trial i exists until duration D_i.
An available-observation mean is

mu_hat(t) =
    sum_i I(D_i >= t) X_i(t)
    -------------------------
        sum_i I(D_i >= t)

The problem is that its interpretation depends on why the epoch ended.
For a self-paced task, this can become approximately

E[X(t) | RT >= t]

because at late times the fast-response trials are already gone.
So the population represented by the average changes with time.
For a spindle, a similar average at late times describes long-lasting spindles, not necessarily the population of all spindles.
For gait, it increasingly describes long strides.
For genuinely missing observations, an available-case estimator may be justified under assumptions about the observation mechanism.
These are different statistical situations.
That is why I don't think np.nanmean() should become the default just because it returns a rectangular result.
nave is another reason to be careful
Suppose:

t = 0.1 s : 100 epochs contribute
t = 0.8 s : 70 epochs contribute
t = 1.5 s : 15 epochs contribute

There is no single honest value of:

evoked.nave

for the whole trajectory.
That is not only display metadata. nave also interacts with noise covariance / inverse scaling in MNE.
So creating a normal Evoked from an available-case ragged average can propagate an incorrect assumption into source estimation.
For this reason I think average() should continue to raise until we decide what output type and statistical contract we actually want.
TFR should be separated into two problems
A per-segment TFR is comparatively straightforward:

X_i(channel, time_i)
        |
        v
P_i(channel, frequency, time_i)

The frequency axis can be common while the time axis remains ragged.
So I think a future non-averaged ragged TFR representation is plausible.
The hard problem is the reduction across epochs.
For phase-varying processes such as gait, temporal registration may be appropriate, but transformation order matters.
For spectral analysis we should preserve physical frequency:

EEG(t)
  |
  v
TFR(f, t)
  |
  v
register the TFR time axis
  |
  v
TFR(f, phase)

rather than

warp EEG(t)
  |
  v
TFR

because stretching/compressing the raw signal changes its physical frequencies.
So I would treat

compute_tfr(..., average=False)

and

compute_tfr(..., average=True)

as scientifically different problems.
PSD also needs an explicit contract
PSD can be estimated independently from individual segments, but unequal duration creates two issues.
First, spectral resolution differs.
Very short segments cannot provide the same low-frequency information as long segments.
A common frequency evaluation grid does not imply common spectral resolution.
Second, averaging requires a weighting decision.
Event-weighted:

P_event(f) = (1 / N) * sum_i P_i(f)

gives every event the same influence.
Duration-weighted:

P_time(f) =
    sum_i D_i P_i(f)
    ----------------
       sum_i D_i

gives every observed second equal influence.
Both can be scientifically legitimate.
They answer different questions.
MNE should therefore not choose one implicitly.
Covariance has the same weighting problem
If samples are pooled:

C_sample =
    sum_i X_i X_i^T
    ----------------
       sum_i T_i

then longer segments contribute proportionally more observations.
An alternative is:

C_epoch = (1 / N) * sum_i C_i

where every epoch covariance gets equal weight.
These are different estimators.
The right choice depends on whether we want the covariance of:
a randomly selected sample;
a randomly selected event;
a condition-balanced population;
or something else.
This is especially important because covariance propagates into:
whitening;
minimum-norm inverse estimation;
beamformers;
CSP;
SPoC;
Mahalanobis-type metrics;
and other downstream methods.
So covariance should be treated as one of the higher-priority research questions rather than just implementing ragged concatenation.
ICA
ICA itself does not fundamentally require equal-length trials.
Ragged epochs could technically be concatenated:

[X_1 | X_2 | X_3 | ...]

But then a 10-second segment contributes 10 times as many observations as a 1-second segment.
For artifact decomposition that may actually be desirable.
For a condition-balanced scientific decomposition it may not be.
My current preference therefore remains:

Fit ICA on continuous Raw before segmentation whenever possible.
We can consider ragged-epoch ICA later once its weighting semantics are explicit.
Source estimation
This seems easier once the spatial operator already exists.
For a fixed inverse operator W:

s_i(t) = W X_i(t)

there is no mathematical requirement that every segment has the same number of time samples.
So minimum-norm / beamformer application could probably eventually produce variable-length source estimates.
The difficult issue is upstream:
how the noise covariance was estimated;
how data covariance was estimated;
how many observations contributed;
whether conditions/durations were weighted appropriately.
So I would separate:

source operator estimation
from
applying an existing source operator to ragged data.
Dipole fitting
Dipole fitting should probably remain blocked until we have a meaningful common representation.
If there is a defensible common physical-time window:

event -> fixed common interval -> Evoked -> dipole

is fine.
If there is a meaningful phase coordinate:

registered phase -> average topography -> phase-specific dipole

may also be meaningful.
But if neither exists, we should not fabricate an Evoked solely so that fit_dipole() accepts the data.
Single-segment / single-trial dipole analysis is a separate possibility, with its own SNR limitations.
Time-resolved decoding and cluster statistics
These are fundamentally common-coordinate methods.
For example, a sliding estimator assumes that:

time point 100

has the same scientific meaning across observations.
Cluster statistics also require a shared lattice / adjacency across observations.
So ragged input should first undergo an explicit transformation such as:
physical-time intersection;
event locking;
response locking;
landmark registration;
phase normalization;
or fixed-dimensional feature extraction.
Without one of those, I think these methods should raise.
Different event types make reductions even more dangerous
Variable-duration datasets may contain:

A -> B -> C -> A -> B -> D

where each type represents a different process.
So even after solving duration, the user should normally select comparable segment types before reduction.
For example:

segments["movement"]

versus

segments["preparation"]

rather than accidentally averaging all intervals together simply because they live in the same container.
Variable time axes are slightly broader than variable duration
Another API point worth keeping in mind:

epoch A: 0 -> 1 s
epoch B: 1 -> 2 s

Both epochs have the same duration, but they do not have the same time axis.
So the deepest property is not strictly:

variable_duration

but something more like:

has_common_times

or

is_time_aligned

This is also relevant to #5794, where the core issue is different per-trial temporal anchors rather than necessarily different numbers of samples.
I don't think we need to solve this API naming question immediately, but we should keep the distinction visible before freezing more public semantics.
Validation plan
Before making the statistical methods ragged-aware, I think we should validate the design against a small set of scientifically different paradigms rather than one dataset.
A. Repeated phase process
Gait cycles
Tests:
variable duration;
known internal landmarks;
TFR-before-registration;
physical frequency preservation;
phase-domain averaging;
stride duration retained as metadata.
B. Naturally bounded transient
Sleep spindle / seizure-like event
Tests:
onset/offset locking;
true event duration;
no interpretation of post-event time as missing data;
optional normalized event progression only when explicitly requested.
C. Response-terminated cognition
Self-paced arithmetic / decision task
Tests:
stimulus-locked representation;
response-locked representation;
demonstrate changing trial population in naive available-time averaging;
model RT explicitly.
This maps directly onto #5612.
D. Genuine partial observation
Synthetic or real truncated recordings
Tests:
underlying complete trajectory known in simulation;
MCAR vs informative truncation;
compare available-case and partial-functional estimators;
explicitly demonstrate when nanmean is valid/invalid.
E. Dense overlapping events
Fixations / words / tone sequences
Tests:
overlapping responses;
repeated use of the same Raw samples;
compare epoch-based representation with regression/deconvolution.
This is related to #5794 and #11480.
F. State sequence
Sleep-stage bouts / anesthesia states
Tests:
variable segment type;
variable dwell duration;
variable number/order of states;
distinguish per-state signal analysis from state-transition modelling.
G. Variable stimulus content
Sentence / speech task
Tests:
variable stimulus duration;
multiple internal landmarks such as words/phonemes;
compare naive 0-100% normalization against landmark or continuous encoding models.
This maps back to #3533.
Checklist for every future downstream method
Before implementing a ragged version, I think we should answer:
What is the estimand?
Does the operation require a shared time coordinate?
Does it require a shared frequency coordinate/resolution?
Does it weight samples or epochs?
Does unequal duration change that weighting?
Is termination natural, behavioral, or missing-data censoring?
Is missingness assumed MCAR/MAR/informative?
Is physical time preserved?
Is physical frequency preserved?
Are internal landmarks available?
Could epochs overlap and therefore duplicate Raw samples?
Does the method depend on covariance?
Does it affect inverse/source estimation?
What happens when contributor count changes with time?
What are the known failure cases?
Does the method reduce exactly to current MNE behaviour when all epochs share the same time axis?
I think this checklist should become our basic criterion for accepting later variable-domain method implementations.
Literature starting points
These seem particularly useful for the statistical contracts rather than merely the container.
Kraus, 2019 — "Inferential procedures for partially observed functional data." Journal of Multivariate Analysis.
DOI: 10.1016/j.jmva.2019.05.002
Relevant for means, covariance, PCA and inference when trajectories are genuinely only partially observed.
Liebl & Rameseder, 2019 — "Partially observed functional data: The case of systematically missing parts."
DOI: 10.1016/j.csda.2018.08.011
Important because available-case estimation can fail when observation duration depends on the underlying trajectory.
Gellar et al., 2014 — "Variable-Domain Functional Regression for Modeling ICU Data."
DOI: 10.1080/01621459.2014.940044
Relevant because it treats domain length itself as meaningful instead of automatically mapping every trajectory onto [0, 1].
Gwin et al., 2011 — "Electrocortical activity is coupled to gait cycle phase during treadmill walking."
DOI: 10.1016/j.neuroimage.2010.08.066
Representative of the gait-phase use case.
Ehinger & Dimigen, 2019 — "Unfold: an integrated toolbox for overlap correction, non-linear modeling, and regression-based EEG analysis."
DOI: 10.7717/peerj.7838
Relevant for overlapping events where creating independent epochs can itself be the wrong model.
Borst & Anderson, 2015 — "The discovery of processing stages: Analyzing EEG data with hidden semi-Markov models."
DOI: 10.1016/j.neuroimage.2014.12.029
Relevant for variable-duration cognitive/state sequences where state duration should be explicitly modelled.
I think we should add more references as we tackle each method family rather than trying to settle all of the statistics in this PR.
How the historical MNE issues fit into this
The previous requests also make more sense when separated this way:
#3533 — variable-duration representation, especially sentences.
#5612 — variable-duration self-paced cognitive intervals and TFR.
#5794 — multiple temporal anchors / per-trial realignment.
#11480 — baseline relative to an event other than the epoch-locking event.
#12315 — PSD from annotation-defined ragged intervals.
#14206 — umbrella issue where we can separate representation from downstream statistical policy.
So I think the old issues are related because they expose missing infrastructure, but they are not asking for one universal estimator.
What I suggest for #14210
For this PR specifically, I think we should keep the scope intentionally conservative.
Include:
explicit variable-duration extraction / representation;
real samples only;
per-epoch bounds and durations;
per-epoch data/time access;
basic indexing and channel operations;
correct browsing;
an explicit rectangular conversion when the user asks for it;
contributor-count information for that conversion;
unchanged behaviour for normal fixed-duration Epochs;
tests for very short/long combinations, dropped epochs, indexing, browsing, etc.
Continue to raise for methods where we have not defined the scientific contract.
Then work through the other method families gradually in #14206 and separate follow-up PRs.
My preference would be roughly:

PR 1
representation + browsing
        |
        v
PR 2
safe per-segment transforms
        |
        v
PR 3
ragged PSD semantics
        |
        v
PR 4
non-averaged ragged TFR
        |
        v
research / validation
        |
        +--> averaging / Evoked
        +--> covariance
        +--> ICA / CSP
        +--> source-estimation contracts
        +--> statistics / decoding

The exact order can change, but I think small method-family PRs with explicit validation are much safer than trying to make the whole Epochs API accept ragged data at once.
So the goal I would use going forward is not:

every MNE function should accept variable-duration epochs.
It is:
Every MNE function that accepts variable-domain data should either have a clearly defined statistical meaning with visible assumptions, or explicitly refuse until the user supplies the missing scientific decision.
That seems to me the safest way to make this genuinely useful without creating behaviour that looks convenient but quietly changes the science.

larsoner reviewed Aug 25, 2026

View reviewed changes

larsoner added this to 2026 MNE-Python Maintainers Sprint Aug 25, 2026

github-project-automation Bot moved this to Ideas in 2026 MNE-Python Maintainers Sprint Aug 25, 2026

larsoner moved this from Ideas to In Progress in 2026 MNE-Python Maintainers Sprint Aug 25, 2026

snesmaeili mentioned this pull request Aug 27, 2026

MNE MoBI mne-tools/mne-denoise#81

Open

snesmaeili force-pushed the epochs-variable-duration branch from e9808c4 to 0d3b433 Compare August 27, 2026 10:21

snesmaeili mentioned this pull request Aug 27, 2026

ENH: support variable-duration epochs in browser mne-tools/mne-qt-browser#452

Draft

snesmaeili added 4 commits August 27, 2026 13:44

snesmaeili force-pushed the epochs-variable-duration branch from 0d3b433 to bed56d5 Compare August 27, 2026 12:01

snesmaeili marked this pull request as ready for review August 27, 2026 12:35

snesmaeili requested review from agramfort, drammock and mscheltienne as code owners August 27, 2026 12:35

larsoner reviewed Aug 27, 2026

View reviewed changes

snesmaeili added 5 commits August 27, 2026 22:17

FIX: skip the metadata crop test when pandas is missing

75d6e5a

`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.

Navigation Menu

Uh oh!

Uh oh!

ENH: allow Epochs to hold trials of different duration - #14210

ENH: allow Epochs to hold trials of different duration#14210
snesmaeili wants to merge 11 commits into
mne-tools:mainfrom
snesmaeili:epochs-variable-duration

snesmaeili commented Aug 25, 2026 •

edited

Loading

Uh oh!

larsoner left a comment

Uh oh!

larsoner Aug 25, 2026

Uh oh!

larsoner Aug 25, 2026

Uh oh!

snesmaeili Aug 25, 2026

Uh oh!

larsoner Aug 25, 2026

Uh oh!

larsoner Aug 25, 2026

Uh oh!

drammock commented Aug 25, 2026

Uh oh!

larsoner left a comment

Uh oh!

larsoner Aug 27, 2026

Uh oh!

snesmaeili commented Aug 28, 2026 •

edited

Loading

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

	tmin: "float \| np.ndarray" = -0.2,
	tmin: float \| np.ndarray = -0.2,



		def _check_variable_bounds(tmin, tmax, n_events):
		"""Normalize ``tmin``/``tmax``, which may be given per event.

Uh oh!

Uh oh!

Conversation

snesmaeili commented Aug 25, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

larsoner left a comment

Choose a reason for hiding this comment

Uh oh!

larsoner Aug 25, 2026

Choose a reason for hiding this comment

Uh oh!

larsoner Aug 25, 2026

Choose a reason for hiding this comment

Uh oh!

snesmaeili Aug 25, 2026

Choose a reason for hiding this comment

Uh oh!

larsoner Aug 25, 2026

Choose a reason for hiding this comment

Uh oh!

larsoner Aug 25, 2026

Choose a reason for hiding this comment

Uh oh!

drammock commented Aug 25, 2026

Uh oh!

larsoner left a comment

Choose a reason for hiding this comment

Uh oh!

larsoner Aug 27, 2026

Choose a reason for hiding this comment

Uh oh!

snesmaeili commented Aug 28, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

snesmaeili commented Aug 25, 2026 •

edited

Loading

snesmaeili commented Aug 28, 2026 •

edited

Loading