[7f9fb8]: / tutorials / evoked / 10_evoked_overview.py

Download this file

346 lines (305 with data), 15.5 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
"""
.. _tut-evoked-class:
===============================================
The Evoked data structure: evoked/averaged data
===============================================
This tutorial covers the basics of creating and working with :term:`evoked` data. It
introduces the :class:`~mne.Evoked` data structure in detail, including how to load,
query, subset, export, and plot data from an :class:`~mne.Evoked` object. For details
on creating an :class:`~mne.Evoked` object from (possibly simulated) data in a
:class:`NumPy array <numpy.ndarray>`, see :ref:`tut-creating-data-structures`.
As usual, we start by importing the modules we need:
"""
# Authors: The MNE-Python contributors.
# License: BSD-3-Clause
# Copyright the MNE-Python contributors.
# %%
import mne
# %%
# Creating ``Evoked`` objects from ``Epochs``
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# :class:`~mne.Evoked` objects typically store EEG or MEG signals that have
# been *averaged* over multiple :term:`epochs`, which is a common technique for
# estimating stimulus-evoked activity. The data in an :class:`~mne.Evoked`
# object are stored in an :class:`array <numpy.ndarray>` of shape
# ``(n_channels, n_times)`` (in contrast to an :class:`~mne.Epochs` object,
# which stores data of shape ``(n_epochs, n_channels, n_times)``). Thus, to
# create an :class:`~mne.Evoked` object, we'll start by epoching some raw data,
# and then averaging together all the epochs from one condition:
root = mne.datasets.sample.data_path() / "MEG" / "sample"
raw_file = root / "sample_audvis_raw.fif"
raw = mne.io.read_raw_fif(raw_file, verbose=False)
events = mne.find_events(raw, stim_channel="STI 014")
# we'll skip the "face" and "buttonpress" conditions to save memory
event_dict = {
"auditory/left": 1,
"auditory/right": 2,
"visual/left": 3,
"visual/right": 4,
}
epochs = mne.Epochs(raw, events, tmin=-0.3, tmax=0.7, event_id=event_dict, preload=True)
evoked = epochs["auditory/left"].average()
del raw # reduce memory usage
# %%
# You may have noticed that MNE informed us that "baseline correction" has been
# applied. This happened automatically during creation of the
# :class:`~mne.Epochs` object, but may also be initiated (or disabled)
# manually. We will discuss this in more detail later.
#
# The information about the baseline period of :class:`~mne.Epochs` is
# transferred to derived :class:`~mne.Evoked` objects to maintain provenance as
# you process your data:
print(f"Epochs baseline: {epochs.baseline}")
print(f"Evoked baseline: {evoked.baseline}")
# %%
# Basic visualization of ``Evoked`` objects
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# We can visualize the average evoked response for left-auditory stimuli using
# the :meth:`~mne.Evoked.plot` method, which yields a butterfly plot of each
# channel type:
evoked.plot()
# %%
# Like the ``plot()`` methods for :meth:`Raw <mne.io.Raw.plot>` and
# :meth:`Epochs <mne.Epochs.plot>` objects,
# :meth:`evoked.plot() <mne.Evoked.plot>` has many parameters for customizing
# the plot output, such as color-coding channel traces by scalp location, or
# plotting the :term:`global field power` alongside the channel traces.
# See :ref:`tut-visualize-evoked` for more information on visualizing
# :class:`~mne.Evoked` objects.
#
#
# Subsetting ``Evoked`` data
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# .. admonition:: Evokeds are not memory-mapped
# :class: sidebar note
#
# :class:`~mne.Evoked` objects use a :attr:`~mne.Evoked.data` *attribute*
# rather than a :meth:`~mne.Epochs.get_data` *method*; this reflects the fact
# that the data in :class:`~mne.Evoked` objects are always loaded into
# memory and never `memory-mapped`_ from their location on disk (because they
# are typically *much* smaller than :class:`~mne.io.Raw` or
# :class:`~mne.Epochs` objects).
#
#
# Unlike :class:`~mne.io.Raw` and :class:`~mne.Epochs` objects,
# :class:`~mne.Evoked` objects do not support selection by square-bracket
# indexing. Instead, data can be subsetted by indexing the
# :attr:`~mne.Evoked.data` attribute:
print(evoked.data[:2, :3]) # first 2 channels, first 3 timepoints
# %%
# To select based on time in seconds, the :meth:`~mne.Evoked.time_as_index`
# method can be useful, although beware that depending on the sampling
# frequency, the number of samples in a span of given duration may not always
# be the same (see the :ref:`time-as-index` section of the tutorial on
# :ref:`Raw data <tut-raw-class>` for details).
#
#
# Selecting, dropping, and reordering channels
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#
# By default, when creating :class:`~mne.Evoked` data from an
# :class:`~mne.Epochs` object, only the primary data channels will be retained:
# ``eog``, ``ecg``, ``stim``, and ``misc`` channel types will be dropped. You
# can control which channel types are retained via the ``picks`` parameter of
# :meth:`epochs.average() <mne.Epochs.average>`, by passing ``'all'`` to
# retain all channels, or by passing a list of integers, channel names, or
# channel types. See the documentation of :meth:`~mne.Epochs.average` for
# details.
#
# If you've already created the :class:`~mne.Evoked` object, you can use the
# :meth:`~mne.Evoked.pick`, :meth:`~mne.Evoked.pick_channels`,
# :meth:`~mne.Evoked.pick_types`, and :meth:`~mne.Evoked.drop_channels` methods
# to modify which channels are included in an :class:`~mne.Evoked` object.
# You can also use :meth:`~mne.Evoked.reorder_channels` for this purpose; any
# channel names not provided to :meth:`~mne.Evoked.reorder_channels` will be
# dropped. Note that *channel* selection methods modify the object in-place, so
# in interactive/exploratory sessions you may want to create a
# :meth:`~mne.Evoked.copy` first.
evoked_eeg = evoked.copy().pick(picks="eeg")
print(evoked_eeg.ch_names)
new_order = ["EEG 002", "MEG 2521", "EEG 003"]
evoked_subset = evoked.copy().reorder_channels(new_order)
print(evoked_subset.ch_names)
# %%
# Similarities among the core data structures
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# :class:`~mne.Evoked` objects have many similarities with :class:`~mne.io.Raw`
# and :class:`~mne.Epochs` objects, including:
#
# - They can be loaded from and saved to disk in ``.fif`` format, and their
# data can be exported to a :class:`NumPy array <numpy.ndarray>` (but through
# the :attr:`~mne.Evoked.data` attribute instead of a ``get_data()``
# method). :class:`Pandas DataFrame <pandas.DataFrame>` export is also
# available through the :meth:`~mne.Evoked.to_data_frame` method.
#
# - You can change the name or type of a channel using
# :meth:`evoked.rename_channels() <mne.Evoked.rename_channels>` or
# :meth:`evoked.set_channel_types() <mne.Evoked.set_channel_types>`.
# Both methods take :class:`dictionaries <dict>` where the keys are existing
# channel names, and the values are the new name (or type) for that channel.
# Existing channels that are not in the dictionary will be unchanged.
#
# - :term:`SSP projector <projector>` manipulation is possible through
# :meth:`~mne.Evoked.add_proj`, :meth:`~mne.Evoked.del_proj`, and
# :meth:`~mne.Evoked.plot_projs_topomap` methods, and the
# :attr:`~mne.Evoked.proj` attribute. See :ref:`tut-artifact-ssp` for more
# information on SSP.
#
# - Like :class:`~mne.io.Raw` and :class:`~mne.Epochs` objects,
# :class:`~mne.Evoked` objects have :meth:`~mne.Evoked.copy`,
# :meth:`~mne.Evoked.crop`, :meth:`~mne.Evoked.time_as_index`,
# :meth:`~mne.Evoked.filter`, and :meth:`~mne.Evoked.resample` methods.
#
# - Like :class:`~mne.io.Raw` and :class:`~mne.Epochs` objects,
# :class:`~mne.Evoked` objects have ``evoked.times``,
# :attr:`evoked.ch_names <mne.Evoked.ch_names>`, and :class:`info <mne.Info>`
# attributes.
#
#
# .. _tut-section-load-evk:
#
# Loading and saving ``Evoked`` data
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# Single :class:`~mne.Evoked` objects can be saved to disk with the
# :meth:`evoked.save() <mne.Evoked.save>` method. One difference between
# :class:`~mne.Evoked` objects and the other data structures is that multiple
# :class:`~mne.Evoked` objects can be saved into a single ``.fif`` file, using
# :func:`mne.write_evokeds`. The :ref:`example data <sample-dataset>`
# includes such a ``.fif`` file: the data have already been epoched and
# averaged, and the file contains separate :class:`~mne.Evoked` objects for
# each experimental condition:
evk_file = root / "sample_audvis-ave.fif"
evokeds_list = mne.read_evokeds(evk_file, verbose=False)
print(evokeds_list)
print(type(evokeds_list))
# %%
# Notice that :func:`mne.read_evokeds` returned a :class:`list` of
# :class:`~mne.Evoked` objects, and each one has an ``evoked.comment``
# attribute describing the experimental condition that was averaged to
# generate the estimate:
for evok in evokeds_list:
print(evok.comment)
# %%
# If you want to load only some of the conditions present in a ``.fif`` file,
# :func:`~mne.read_evokeds` has a ``condition`` parameter, which takes either a
# string (matched against the comment attribute of the evoked objects on disk),
# or an integer selecting the :class:`~mne.Evoked` object based on the order
# it is stored in the file. Passing lists of integers or strings is also
# possible. If only one object is selected, the :class:`~mne.Evoked` object
# will be returned directly (rather than inside a list of length one):
right_vis = mne.read_evokeds(evk_file, condition="Right visual")
print(right_vis)
print(type(right_vis))
# %%
# Previously, when we created an :class:`~mne.Evoked` object by averaging
# epochs, baseline correction was applied by default when we extracted epochs
# from the `~mne.io.Raw` object (the default baseline period is ``(None, 0)``,
# which ensures zero mean for times before the stimulus event). In contrast, if
# we plot the first :class:`~mne.Evoked` object in the list that was loaded
# from disk, we'll see that the data have not been baseline-corrected:
evokeds_list[0].plot(picks="eeg")
# %%
# This can be remedied by either passing a ``baseline`` parameter to
# :func:`mne.read_evokeds`, or by applying baseline correction after loading,
# as shown here:
# Original baseline (none set)
print(f"Baseline after loading: {evokeds_list[0].baseline}")
# Apply a custom baseline correction
evokeds_list[0].apply_baseline((None, 0))
print(f"Baseline after calling apply_baseline(): {evokeds_list[0].baseline}")
# Visualize the evoked response
evokeds_list[0].plot(picks="eeg")
# %%
# Notice that :meth:`~mne.Evoked.apply_baseline` operated in-place. Similarly,
# :class:`~mne.Evoked` objects may have been saved to disk with or without
# :term:`projectors <projector>` applied; you can pass ``proj=True`` to the
# :func:`~mne.read_evokeds` function, or use the :meth:`~mne.Evoked.apply_proj`
# method after loading.
#
#
# Combining ``Evoked`` objects
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# One way to pool data across multiple conditions when estimating evoked
# responses is to do so *prior to averaging* (recall that MNE-Python can select
# based on partial matching of epoch labels separated by ``/``; see
# :ref:`tut-section-subselect-epochs` for more information):
left_right_aud = epochs["auditory"].average()
left_right_aud
# %%
# This approach will weight each epoch equally and create a single
# :class:`~mne.Evoked` object. Notice that the printed representation includes
# ``(average, N=145)``, indicating that the :class:`~mne.Evoked` object was
# created by averaging across 145 epochs. In this case, the event types were
# fairly close in number:
left_aud = epochs["auditory/left"].average()
right_aud = epochs["auditory/right"].average()
print([evok.nave for evok in (left_aud, right_aud)])
# %%
# However, this may not always be the case. If for statistical reasons it is
# important to average *the same number* of epochs from different conditions,
# you can use :meth:`~mne.Epochs.equalize_event_counts` prior to averaging.
#
# Another approach to pooling across conditions is to create separate
# :class:`~mne.Evoked` objects for each condition, and combine them afterwards.
# This can be accomplished with the function :func:`mne.combine_evoked`, which
# computes a weighted sum of the :class:`~mne.Evoked` objects given to it. The
# weights can be manually specified as a list or array of float values, or can
# be specified using the keyword ``'equal'`` (weight each :class:`~mne.Evoked`
# object by :math:`\frac{1}{N}`, where :math:`N` is the number of
# :class:`~mne.Evoked` objects given) or the keyword ``'nave'`` (weight each
# :class:`~mne.Evoked` object proportional to the number of epochs averaged
# together to create it):
left_right_aud = mne.combine_evoked([left_aud, right_aud], weights="nave")
assert left_right_aud.nave == left_aud.nave + right_aud.nave
# %%
# Note that the ``nave`` attribute of the resulting :class:`~mne.Evoked` object
# will reflect the *effective* number of averages, and depends on both the
# ``nave`` attributes of the contributing :class:`~mne.Evoked` objects and the
# weights with which they are combined. Keeping track of effective ``nave`` is
# important for inverse imaging, because ``nave`` is used to scale the noise
# covariance estimate, which in turn affects the magnitude of estimated source
# activity (see :ref:`minimum_norm_estimates` for more information, especially
# the :ref:`whitening_and_scaling` section). Note that
# :func:`mne.grand_average` does *not* adjust ``nave`` to reflect the effective
# number of averaged epochs; it simply sets ``nave`` to the number of *evokeds*
# that were averaged together. For this reason, it is best to use
# :func:`mne.combine_evoked` rather than :func:`mne.grand_average` if you
# intend to perform inverse imaging on the resulting :class:`~mne.Evoked`
# object.
#
#
# Other uses of ``Evoked`` objects
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#
# Although the most common use of :class:`~mne.Evoked` objects is to store
# *averages* of epoched data, there are a few other uses worth noting here.
# First, the method :meth:`epochs.standard_error() <mne.Epochs.standard_error>`
# will create an :class:`~mne.Evoked` object (just like
# :meth:`epochs.average() <mne.Epochs.average>` does), but the data in the
# :class:`~mne.Evoked` object will be the standard error across epochs instead
# of the average. To indicate this difference, :class:`~mne.Evoked` objects
# have a :attr:`~mne.Evoked.kind` attribute that takes values ``'average'`` or
# ``'standard error'`` as appropriate.
#
# Another use of :class:`~mne.Evoked` objects is to represent *a single trial
# or epoch* of data, usually when looping through epochs. This can be easily
# accomplished with the :meth:`epochs.iter_evoked() <mne.Epochs.iter_evoked>`
# method, and can be useful for applications where you want to do something
# that is only possible for :class:`~mne.Evoked` objects. For example, here
# we use the :meth:`~mne.Evoked.get_peak` method (which is not available for
# :class:`~mne.Epochs` objects) to get the peak response in each trial:
for ix, trial in enumerate(epochs[:3].iter_evoked()):
channel, latency, value = trial.get_peak(ch_type="eeg", return_amplitude=True)
latency = int(round(latency * 1e3)) # convert to milliseconds
value = int(round(value * 1e6)) # convert to µV
print(f"Trial {ix}: peak of {value} µV at {latency} ms in channel {channel}")
# %%
# .. REFERENCES
#
# .. _`memory-mapped`: https://en.wikipedia.org/wiki/Memory-mapped_file