Question: decim after reading epochs

[teonbrooks]

Is there a way to decimate Epochs read in with mne.read_epochs? There's the resampling method but it's doing something different than mne.Epochs(..., decim)? If there isn't, should we add the parameter to mne.read_epochs?

[larsoner]

I think we should just add a .decimate(n) method

[agramfort]

Why not decimating when creating epochs?

[teonbrooks]

I forgot :( and macbook is weak :((

[agramfort]

Rerunning will be faster than doing a PR ;)

[teonbrooks]

you haven't met my 5 yo two-core system ;) but aside from my data, do you think that would be a useful feature? for instance, if you have data generated with EpochsArray from some mat file, or read_epochs*, those files can't be decimated post-creation

[agramfort]

I am -0.5 on this. There should be ideally 1 way to do things. We offer this at Epochs creation.

my 2c

[jona-sassenhagen]

Seems like a useful feature for example if you want to 1. have many sampling points for a data-hungry preprocessing stage like ICA, but want to 2. run a later analysis step on decimated epochs to save time.

Dropping epochs from an Epochs object by amplitude is a similar thing.

Maybe Epochs could take not only Raw, but also Epochs instances, so the init stuff could simply be reused?

[dengemann]

@agramfort , @fraimondo and I but also @kingjr are using privately functions to decimate after epochs creations, there are sevaral use cases.

I'd be +50 for adding a tested epochs.decim function.

Fede and I are using it for experiments where we compute down-sampled epochs (that we store) and later apply several decimation factors.

2015-06-02 21:46 GMT+02:00 jona-sassenhagen notifications@github.com:

Seems like a useful feature for example if you want to 1. have many sampling points for a data-hungry preprocessing stage like ICA, but want to

2. run a later analysis step on decimated epochs to save time.

Dropping epochs from an Epochs object by amplitude is a similar thing.

Maybe Epochs could take not only Raw, but also Epochs instances, so the init stuff could simply be reused?

— Reply to this email directly or view it on GitHub https://github.com/mne-tools/mne-python/issues/2178#issuecomment-108075037 .

[agramfort]

ok you win...

[kingjr]

Yes, I'm also super +1 for epoch decimation. The most classic scenario is people sending me heavy epoched data.

[teonbrooks]

anybody want to share their private method and/or start a pr?

[kingjr]

I just did this:

def decimate(inst, decim):
    """Fast MNE object decimation"""
    from mne.io.base import _BaseRaw
    from mne.epochs import _BaseEpochs
    if isinstance(inst, _BaseRaw):
        inst._data = inst._data[:, ::decim]
        inst.info['sfreq'] /= decim
        inst._first_samps /= decim
        inst._last_samps /= decim
    elif isinstance(inst, _BaseEpochs):
        inst._data = inst._data[:, :, ::decim]
        inst.info['sfreq'] /= decim
        inst.times = inst.times[::decim]
    return inst

[fraimondo]

We did something like this, although it's part of a bigger module that allows us to subsample in spatial or temporal manners.

def subsample(inst, n_channels=None, decim=None):
    _subsample_channels(inst, n_channels)
    if isinstance(inst, mne.io.base._BaseRaw):
        _subsample_raw(inst, decim=decim)
    elif isinstance(inst, mne.epochs._BaseEpochs):
        _subsample_epochs(inst, decim=decim)

def _subsample_epochs(epochs, n_channels=None, decim=None):
    if decim is not None:
        if decim > 1:  # Do not decimate by less than 1
            logger.info('Decimating by {}'.format(decim))
            epochs.times = epochs.times[::decim]
            epochs._data = epochs._data[..., ::decim]
            epochs.info['sfreq'] /= decim

def _subsample_raw(raw, n_channels=None, decim=None):
    if decim is not None:
        if decim > 1:  # Do not decimate by less than 1
            logger.info('Decimating by {}'.format(decim))
            events = mne.find_events(raw)
            events[:, 0] = np.ceil(events[:, 0] / float(decim)).astype(int)
            pick_stim = mne.pick_types(raw.info, meg=False, eeg=False, stim=True)
            raw._last_samps /= decim
            raw._data = raw._data[:, ::decim]
            raw._data[pick_stim].fill(0.0)
            raw._data[pick_stim, events[:, 0]] = events[:, 2]
            assert np.all(raw._data[pick_stim, events[:, 0]] == events[:, 2])
            raw.info['sfreq'] /= decim

[jona-sassenhagen]

Again, wouldn't setting up mne.Epochs to allow an Epochs instance solve a bunch of problems at once? You could also re-baseline, re-reject, crop ... all without introducing any new public API. It would all be internal.

[dengemann]

Again, wouldn't setting up mne.Epochs to allow an Epochs instance solve a bunch of problems at once? You could also re-baseline, re-reject, crop ... all without introducing any new public API. It would all be internal.

@jona-sassenhagen could you rephrase this? I'm not sure I got your proposal. You might suggest that we have an ResampleMixin that's inhereted by Raw, Epochs et al. with a decimate method? Cropping is already possible everywhere via a method. Re-baselining would indeed be cool, but it would require storing, for example the mean-vector, for a mean baseline.

[jona-sassenhagen]

Pretty sure I got the idea from a code snippet by you. Basically, imagine the API for rebaselining to be this: new_epochs = mne.Epochs(old_epochs, baseline=(-3, 0))

and for decimating: new_epochs = mne.Epochs(old_epochs, decim=2)

and for rejecting a second time: new_epochs = mne.Epochs(old_epochs, reject=dict("eeg":1000))

This way, no new keywords, functions, ... are introduced.

[dengemann]

new_epochs = mne.Epochs(old_epochs, baseline=(-3, 0))

will not work if already baselined

new_epochs = mne.Epochs(old_epochs, decim=2)

this might work.

new_epochs = mne.Epochs(old_epochs, reject=dict("eeg":1000))

we have an open issue for a method. EpochsArray can do it.

But let me think about what you propose for a moment.

[agramfort]

new_epochs = mne.Epochs(old_epochs, baseline=(-3, 0))

will not work if already baselined

why not? I think it would work.

My only fear is that the Epochs class is already a beast and relunctant to make it even more complex.

[teonbrooks]

maybe it's good time to refactor the Epochs like @Eric89GXL suggested before

[dengemann]

maybe it's good time to refactor the Epochs like @Eric89GXL suggested before

I'm -0.5 on big refactories. They create lot's of work and consume time and the gain is arguable at this stage. You need to convince me.

[larsoner]

The epochs class hierarchy is a mess. Take a look at which classes have which methods and what inherits from what (BaseEpochs, Epochs, EpochsArray, etc.) and you'll probably be convinced :) So I am +1 on refactoring. It should lead to simpler code that is easier to maintain.

[larsoner]

And IIRC this bad hierarchy has led to some recent bugs, I assume there are others hiding in there

[agramfort]

Eric are you up for it?

If so before starting can you sketch the plan?

[larsoner]

The plan would be something like:

1. move all of the methods from Epochs to BaseEpochs except for loading (i.e., Epochs means EpochsFIF following the raw analogy)
2. make EpochsArray inherit from BaseEpochs instead of Epochs
3. make EpochsKit inherit from BaseEpochs instead of EpochsArray
4. make sure RtEpochs works (@mainakjas's can you help?) since it won't support all functions. Although maybe this isn't really an issue because we can just set the preloaded parameter False and that will get rid of some of the ones that won't work like resample or savgol.

[larsoner]

So I guess the only part in my mind that isn't clear is what to do with RtEpochs, since it only supports a subset of the methods that all other Epochs instances support. Should I just override methods with NotImplementedError, or do we want BaseEpochs that RtEpochs inherts from, and BaseEpochsFull or so that inherts from BaseEpochs, that all other Epochs classes inherit from (and BaseEpochsFull contains all methods that RtEpochs can't handle)?

[larsoner]

@agramfort this also has a rough plan in it, @dengemann it has some justifications:

https://github.com/mne-tools/mne-python/issues/1982

See also these discussions @dengemann:

https://github.com/mne-tools/mne-python/issues/1633#issuecomment-91924011 https://github.com/mne-tools/mne-python/pull/1590#discussion_r28183638

[agramfort]

or have a light BaseEpochs and an EpochsMixin class that implements all the methods not available for RtEpochs and have Epochs inherit from BaseEpochs and EpochsMixin

[larsoner]

@agramfort in that case Epochs, EpochsArray, and EpochsKit will all inherit from both BaseEpochs and EpochsMixin, okay?

[jona-sassenhagen]

So...are you planning on implementing e.g. Post-hoc decimation in the way I had proposed?

On Sat, Jun 6, 2015 at 6:41 AM -0700, "Eric Larson" notifications@github.com wrote:

@agramfort in that case Epochs, EpochsArray, and EpochsKit will all inherit from both BaseEpochs and EpochsMixin, okay?

— Reply to this email directly or view it on GitHub.

[larsoner]

I'm fine with it, I think @agramfort was reluctant to have it

[teonbrooks]

@jona-sassenhagen I think it would be a method and not just modifying the Epochs class to take an Epochs instance.

+1 for the light _BaseEpochs that all the epochs subclasses use and the EpochsMixin for the additional methods

[mainakjas]

I can help with RtEpochs definitely but I think Martin knows it the best. Should be a good exercise for me :)

[dengemann]

Ok for cleaning up internals. I'm not such a big fan of Epochs that takes an epochs object in it's constructor. A re_epoch function would be cleaner to me, if needed (could then call all methods at once). I would personally be happy about a simple method. It will also have less overhead than calling the constructor.

[larsoner]

+1 for method, I'll try to do it during refactoring