-
Notifications
You must be signed in to change notification settings - Fork 1.3k
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
ENH: dict features for Epochs objects #229
Changes from all commits
95a7f75
4ad8e3f
14d32ed
2f2483a
2755b1c
d0db669
7818662
61e6fbd
9b8e023
d85edcc
4e8f22f
dff3fdc
c1fbf31
f17b4a4
f6cefe6
3381f0c
451d701
a6fde9b
05790d6
b60033e
609728d
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,69 @@ | ||
""" | ||
================================================ | ||
Compare Evoked Reponses for Different Conditions | ||
================================================ | ||
|
||
In this example, an Epochs object for visual and | ||
auditory responses is created. Both conditions | ||
are then accessed by their respective names to | ||
create a sensor layout plot of the related | ||
evoked responses. | ||
|
||
""" | ||
|
||
# Authors: Denis Engemann <d.engemann@fz-juelich.de> | ||
# Alexandre Gramfort <gramfort@nmr.mgh.harvard.edu> | ||
# | ||
# License: BSD (3-clause) | ||
|
||
print __doc__ | ||
|
||
import pylab as pl | ||
import mne | ||
|
||
from mne.fiff import Raw, pick_types | ||
from mne.layouts import read_layout | ||
from mne.viz import plot_topo | ||
from mne.datasets import sample | ||
data_path = sample.data_path('.') | ||
|
||
############################################################################### | ||
# Set parameters | ||
raw_fname = data_path + '/MEG/sample/sample_audvis_filt-0-40_raw.fif' | ||
event_fname = data_path + '/MEG/sample/sample_audvis_filt-0-40_raw-eve.fif' | ||
event_id = 1 | ||
tmin = -0.2 | ||
tmax = 0.5 | ||
|
||
# Setup for reading the raw data | ||
raw = Raw(raw_fname) | ||
events = mne.read_events(event_fname) | ||
|
||
# Set up pick list: MEG + STI 014 - bad channels (modify to your needs) | ||
include = [] # or stim channels ['STI 014'] | ||
exclude = raw.info['bads'] # bads | ||
|
||
# Set up amplitude-peak rejection values for MEG channels | ||
reject = dict(grad=4000e-13, mag=4e-12) | ||
|
||
# pick MEG channels | ||
picks = pick_types(raw.info, meg=True, eeg=False, stim=False, eog=True, | ||
include=include, exclude=exclude) | ||
|
||
# Create epochs including different events | ||
epochs = mne.Epochs(raw, events, dict(audio_l=1, visual_r=3), tmin, tmax, | ||
picks=picks, baseline=(None, 0), reject=reject) | ||
|
||
# Generate list of evoked objects from conditions names | ||
evokeds = [epochs[name].average() for name in 'audio_l', 'visual_r'] | ||
|
||
############################################################################### | ||
# Show topography for two different conditions | ||
|
||
layout = read_layout('Vectorview-all.lout') | ||
|
||
pl.close('all') | ||
pl.figure() | ||
title = 'MNE sample data - left auditory and visual' | ||
plot_topo(evokeds, layout, color=['y', 'g'], title=title) | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. apparently there is a bug if you don't give color then only one evoked is plotted. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I'm addings some cases + guessing. Is there btw. a prefered order? 'w' 'y' 'g' 'r' 'b' ? if not I would just On Fri, Dec 7, 2012 at 8:54 AM, Alexandre Gramfort <notifications@github.com
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more.
itertools.cycle on ['w' 'y' 'g' 'r' 'b'] is good to me. indeed. I can reproduce the swiping… darn... There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. see recent commit On Fri, Dec 7, 2012 at 10:26 AM, Alexandre Gramfort <
|
||
pl.show() |
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -30,6 +30,7 @@ | |
from .filter import resample | ||
from .event import _read_events_fif | ||
from . import verbose | ||
from .fixes import in1d | ||
|
||
|
||
class Epochs(object): | ||
|
@@ -43,8 +44,13 @@ class Epochs(object): | |
events : array, of shape [n_events, 3] | ||
Returned by the read_events function | ||
|
||
event_id : int | None | ||
The id of the event to consider. If None all events are used. | ||
event_id : int | dict | None | ||
The id of the event to consider. If dict, | ||
the keys can later be used to acces associated events. Example: | ||
dict(auditory=1, visual=3). If int, a dict will be created with | ||
the id as string. If None, all events will be used with | ||
and a dict is created with string integer names corresponding | ||
to the event id integers. | ||
|
||
tmin : float | ||
Start time before event | ||
|
@@ -105,6 +111,9 @@ class Epochs(object): | |
info: dict | ||
Measurement info | ||
|
||
event_id : dict | ||
Names of of conditions corresponding to event_ids. | ||
|
||
ch_names: list of string | ||
List of channels' names | ||
|
||
|
@@ -121,7 +130,7 @@ class Epochs(object): | |
get_data() : self | ||
Return all epochs as a 3D array [n_epochs x n_channels x n_times]. | ||
|
||
average(picks=None) : self | ||
average(picks=None) : Evoked | ||
Return Evoked object containing averaged epochs as a | ||
2D array [n_channels x n_times]. | ||
|
||
|
@@ -139,6 +148,14 @@ class Epochs(object): | |
resample() : self, int, int, int, string or list | ||
Resample preloaded data. | ||
|
||
as_data_frame() : DataFrame | ||
Export Epochs object as Pandas DataFrame for subsequent statistical | ||
analyses. | ||
|
||
to_nitime() : TimeSeries | ||
Export Epochs object as nitime TimeSeries object for subsequent | ||
analyses. | ||
|
||
Notes | ||
----- | ||
For indexing and slicing: | ||
|
@@ -148,6 +165,13 @@ class Epochs(object): | |
epochs[idx] : Epochs | ||
Return Epochs object with a subset of epochs (supports single | ||
index and python style slicing) | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Actually it also supports numpy style indexing for the first dimension (events), doesn't it? There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. What do you mean by numpy style? If you mean fancy indexing, yes. On Thu, Dec 6, 2012 at 5:21 PM, Christian Brodbeck <notifications@github.com
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. yes, indexing with an array of ints or booleans |
||
|
||
For subset selection using categorial labels: | ||
|
||
epochs['name'] : Epochs | ||
Return Epochs object with a subset of epochs corresponding to an | ||
experimental condition as specified by event_id. | ||
|
||
""" | ||
@verbose | ||
def __init__(self, raw, events, event_id, tmin, tmax, baseline=(None, 0), | ||
|
@@ -159,10 +183,21 @@ def __init__(self, raw, events, event_id, tmin, tmax, baseline=(None, 0), | |
|
||
self.raw = raw | ||
self.verbose = raw.verbose if verbose is None else verbose | ||
self.event_id = event_id | ||
self.name = name | ||
if isinstance(event_id, dict): | ||
if not all([isinstance(v, int) for v in event_id.values()]): | ||
raise ValueError('Event IDs must be of type integer') | ||
if not all([isinstance(k, str) for k in event_id]): | ||
raise ValueError('Event names must be of type str') | ||
self.event_id = event_id | ||
elif isinstance(event_id, int): | ||
self.event_id = {str(event_id): event_id} | ||
elif event_id is None: | ||
self.event_id = dict((str(e), e) for e in np.unique(events[:, 2])) | ||
else: | ||
raise ValueError('event_id must be dict or int.') | ||
self.tmin = tmin | ||
self.tmax = tmax | ||
self.name = name | ||
self.keep_comp = keep_comp | ||
self.dest_comp = dest_comp | ||
self.baseline = baseline | ||
|
@@ -209,10 +244,9 @@ def __init__(self, raw, events, event_id, tmin, tmax, baseline=(None, 0), | |
|
||
# Select the desired events | ||
self.events = events | ||
if event_id is not None: | ||
selected = np.logical_and(events[:, 1] == 0, | ||
events[:, 2] == event_id) | ||
self.events = self.events[selected] | ||
overlap = in1d(events[:, 2], self.event_id.values()) | ||
selected = np.logical_and(events[:, 1] == 0, overlap) | ||
self.events = events[selected] | ||
|
||
n_events = len(self.events) | ||
|
||
|
@@ -463,26 +497,49 @@ def __repr__(self): | |
def __getitem__(self, key): | ||
"""Return an Epochs object with a subset of epochs | ||
""" | ||
|
||
if self.preload: | ||
data = self._data | ||
del self._data | ||
|
||
epochs = self.copy() | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. If I'm not mistaken this produces a deepcopy of all the data even before indexing, which would not be very efficient. Also, should users modify the Epochs' ._data? If not, why not let numpy decide when to create deep copies of the data (and internally only deepcopy before actually modifying the data)? There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. you're right. But otherwise it means doing the copy manually which is more work to do right. it might work with a trick like:
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. And hope that Eventually it might be cleaner to have an Epochs class that can be initialized with its attributes... and that might also be useful for the result of data transformations such as time-frequency? There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Maybe but it's a huge API break and I think we can live with that now don't you think? There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Yes, hence "eventually" :) If you have a NewEpochs class that initializes with data, then the current Epochs could be transformed into a function that create such objects to keep the api backwards compatible (or subclass to make type-checking work as well). There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more.
self.copy() makes a deep copy of self, including self._data. For consistency with numpy this seems a good idea. However, when indexing this entails a lot of unnecessary copying (especially since you're going to discard most of the copied data anyways when you index). There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Maybe temporarily remove self._data in getitem before calling self.copy()? There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I see your point, but I'm still a bit reluctant to add something like an On Thu, Dec 6, 2012 at 5:37 PM, Christian Brodbeck <notifications@github.com
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Yes sounds reasonable. More than having a private ._copy with index arg... On Thu, Dec 6, 2012 at 5:40 PM, Christian Brodbeck <notifications@github.com
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Added the version you proposed plus fixed a missing assignment-statement. On Thu, Dec 6, 2012 at 5:41 PM, Denis-Alexander Engemann <
|
||
|
||
if self.preload: | ||
self._data, epochs._data = data, data | ||
|
||
if not self._bad_dropped: | ||
warnings.warn("Bad epochs have not been dropped, indexing will be " | ||
"inaccurate. Use drop_bad_epochs() or preload=True") | ||
warnings.warn("Bad epochs have not been dropped, indexing will" | ||
" be inaccurate. Use drop_bad_epochs() or" | ||
" preload=True") | ||
|
||
if isinstance(key, str): | ||
if key not in self.event_id: | ||
raise KeyError('Event "%s" is not in Epochs.' % key) | ||
|
||
epochs = cp.copy(self) # XXX : should use deepcopy but breaks ... | ||
epochs.events = np.atleast_2d(self.events[key]) | ||
key_match = self.events[0:, 2] == self.event_id[key] | ||
epochs.events = self.events[key_match] | ||
|
||
if self.preload: | ||
select = epochs.events[:, 0] | ||
|
||
epochs.name = (key if epochs.name == 'Unknown' | ||
else 'epochs_%s' % key) | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. copying looks good to me now. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Great, thanks heaps for the detailed feedback. On Thu, Dec 6, 2012 at 6:21 PM, Christian Brodbeck <notifications@github.com
|
||
else: | ||
epochs.events = np.atleast_2d(self.events[key]) | ||
if self.preload: | ||
select = key if isinstance(key, slice) else np.atleast_1d(key) | ||
|
||
if self.preload: | ||
if isinstance(key, slice): | ||
epochs._data = self._data[key] | ||
else: | ||
key = np.atleast_1d(key) | ||
epochs._data = self._data[key] | ||
epochs._data = epochs._data[select] | ||
|
||
return epochs | ||
|
||
def average(self, picks=None): | ||
"""Compute average of epochs | ||
|
||
Parameters | ||
---------- | ||
|
||
picks: None | array of int | ||
If None only MEG and EEG channels are kept | ||
otherwise the channels indices in picks are kept. | ||
|
@@ -492,6 +549,7 @@ def average(self, picks=None): | |
evoked : Evoked instance | ||
The averaged epochs | ||
""" | ||
|
||
return self._compute_mean_or_stderr(picks, 'ave') | ||
|
||
def standard_error(self, picks=None): | ||
|
@@ -512,10 +570,8 @@ def standard_error(self, picks=None): | |
|
||
def _compute_mean_or_stderr(self, picks, mode='ave'): | ||
"""Compute the mean or std over epochs and return Evoked""" | ||
if mode == 'stderr': | ||
_do_std = True | ||
else: | ||
_do_std = False | ||
|
||
_do_std = True if mode == 'stderr' else False | ||
evoked = Evoked(None) | ||
evoked.info = cp.deepcopy(self.info) | ||
n_channels = len(self.ch_names) | ||
|
@@ -610,7 +666,7 @@ def crop(self, tmin=None, tmax=None, copy=False): | |
tmask = (self.times >= tmin) & (self.times <= tmax) | ||
tidx = np.where(tmask)[0] | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Do we need the copy? Epochs never modifies the raw object, so this seems unnecessary and inflates memory consumption. Especially when raw has been preloaded and we do e.g. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Ok, a swe wrote the Epochs.copy I wasn't sure whether there might be cases There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I agree maybe we can live without the self.raw.copy() That lead to bugs with raw gets modified between 2 sets of epochs but that can be considered bad practice There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. why not copy_raw argument but it will be unaccessible but the user in the getitem There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. additional copy_raw attribute then for this case? On Thu, Dec 6, 2012 at 3:46 PM, Alexandre Gramfort <notifications@github.com
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. But Epochs doesn't do anything to Raw, it just reads from it. So never copying it should be perfectly fine.. The danger with introducing too many switches like this is that there is an increasing number of combinations that can result in a failure and testing every one of them becomes a combinatorial nightmare. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. old man wisdom :) I agree -1 on copy_raw. It's not necessary as not exposed anyway. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. On Thu, Dec 6, 2012 at 3:47 PM, Martin Luessi notifications@github.comwrote:
Ok, convinced.
Yes, totally agree, before switching to Haskell let's keep up good Wo we will write an epochs.copy that does not copy the raw. Is this what we
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more.
If you want to warn users about to do that, the raw object could have a _has_dependents attribute that is True if an Epochs object has already been created from it, and could issue a warning when modifying the raw data? |
||
|
||
this_epochs = self if not copy else cp.deepcopy(self) | ||
this_epochs = self if not copy else self.copy() | ||
this_epochs.tmin = this_epochs.times[tidx[0]] | ||
this_epochs.tmax = this_epochs.times[tidx[-1]] | ||
this_epochs.times = this_epochs.times[tmask] | ||
|
@@ -643,9 +699,12 @@ def resample(self, sfreq, npad=100, window='boxcar'): | |
def copy(self): | ||
""" Return copy of Epochs instance | ||
""" | ||
raw = self.raw.copy() | ||
raw = self.raw | ||
del self.raw | ||
new = deepcopy(self) | ||
self.raw = raw | ||
new.raw = raw | ||
|
||
return new | ||
|
||
def save(self, fname): | ||
|
@@ -704,7 +763,6 @@ def save(self, fname): | |
end_block(fid, FIFF.FIFFB_MEAS) | ||
end_file(fid) | ||
|
||
|
||
def as_data_frame(self, frame=True): | ||
"""Get the epochs as Pandas panel of data frames | ||
|
||
|
@@ -1011,7 +1069,7 @@ def read_epochs(fname, proj=True, verbose=None): | |
epochs._projector, epochs.info = setup_proj(info) | ||
epochs.ch_names = info['ch_names'] | ||
epochs.baseline = baseline | ||
epochs.event_id = int(np.unique(events[:, 2])) | ||
epochs.event_id = dict((str(e), e) for e in np.unique(events[:, 2])) | ||
fid.close() | ||
|
||
return epochs | ||
|
@@ -1038,7 +1096,7 @@ def bootstrap(epochs, random_state=None): | |
'in the constructor.') | ||
|
||
rng = check_random_state(random_state) | ||
epochs_bootstrap = cp.deepcopy(epochs) | ||
epochs_bootstrap = epochs.copy() | ||
n_events = len(epochs_bootstrap.events) | ||
idx = rng.randint(0, n_events, n_events) | ||
epochs_bootstrap = epochs_bootstrap[idx] | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
nice :)