[med-svn] [python-mne] 311/376: API: compute_covariance now takes as input Epochs
Yaroslav Halchenko
debian at onerussian.com
Fri Nov 27 17:23:14 UTC 2015
This is an automated email from the git hooks/post-receive script.
yoh pushed a commit to annotated tag v0.1
in repository python-mne.
commit 404f762fe621cfad70d1ba37a91e04765f46b962
Author: Alexandre Gramfort <alexandre.gramfort at inria.fr>
Date: Thu Jul 7 10:19:59 2011 +0200
API: compute_covariance now takes as input Epochs
---
mne/__init__.py | 2 +-
mne/cov.py | 104 +++------------------
mne/event.py | 24 +++++
mne/fiff/tests/data/process_raw.sh | 4 +
mne/fiff/tests/data/test-cov.fif | Bin 547234 -> 547234 bytes
mne/fiff/tests/data/test-km-cov.fif | Bin 0 -> 547234 bytes
mne/fiff/tests/data/test.cov | 8 +-
.../tests/data/{test.cov => test_keepmean.cov} | 9 +-
mne/tests/test_cov.py | 40 +++++---
9 files changed, 79 insertions(+), 112 deletions(-)
diff --git a/mne/__init__.py b/mne/__init__.py
index 79fa389..5a01bf5 100644
--- a/mne/__init__.py
+++ b/mne/__init__.py
@@ -2,7 +2,7 @@ __version__ = '0.1.git'
from .cov import read_cov, write_cov, write_cov_file, Covariance, \
compute_raw_data_covariance, compute_covariance
-from .event import read_events, write_events, find_events
+from .event import read_events, write_events, find_events, merge_events
from .forward import read_forward_solution
from .source_estimate import read_stc, write_stc, SourceEstimate, morph_data, \
spatio_temporal_src_connectivity, \
diff --git a/mne/cov.py b/mne/cov.py
index c4f6149..9fd1215 100644
--- a/mne/cov.py
+++ b/mne/cov.py
@@ -20,7 +20,7 @@ from .fiff.write import start_block, end_block, write_int, write_name_list, \
from .fiff.proj import write_proj, make_projector
from .fiff import fiff_open
from .fiff.pick import pick_types, channel_indices_by_type
-from .epochs import Epochs, _is_good
+from .epochs import _is_good
def rank(A, tol=1e-8):
@@ -339,43 +339,6 @@ def read_cov(fid, node, cov_kind):
###############################################################################
# Estimate from data
-
-# def _estimate_compute_covariance_from_epochs(epochs, bmin, bmax, reject, flat,
-# keep_sample_mean):
-# """Estimate noise covariance matrix from epochs
-# """
-# picks_no_eog = pick_types(epochs.info, meg=True, eeg=True, eog=False)
-# n_channels = len(picks_no_eog)
-# ch_names = [epochs.ch_names[k] for k in picks_no_eog]
-# data = np.zeros((n_channels, n_channels))
-# n_samples = 0
-# if bmin is None:
-# bmin = epochs.times[0]
-# if bmax is None:
-# bmax = epochs.times[-1]
-# bmask = (epochs.times >= bmin) & (epochs.times <= bmax)
-#
-# idx_by_type = channel_indices_by_type(epochs.info)
-#
-# for e in epochs:
-#
-# if not _is_good(e, epochs.ch_names, idx_by_type, reject, flat):
-# print "Artefact detected in epoch"
-# continue
-#
-# e = e[picks_no_eog]
-# mu = e[:, bmask].mean(axis=1)
-# e -= mu[:, None]
-# if not keep_sample_mean:
-# e -= np.mean(e, axis=0)
-# data += np.dot(e, e.T)
-# n_samples += e.shape[1]
-#
-# print "Number of samples used : %d" % n_samples
-# print '[done]'
-# return data, n_samples, ch_names
-
-
def compute_raw_data_covariance(raw, tmin=None, tmax=None, tstep=0.2,
reject=None, flat=None, keep_sample_mean=True):
"""Estimate noise covariance matrix from a continuous segment of raw data
@@ -466,51 +429,16 @@ def compute_raw_data_covariance(raw, tmin=None, tmax=None, tstep=0.2,
return cov
-def compute_covariance(raw, events, event_ids, tmin, tmax,
- baseline=(None, None), reject=None, flat=None,
- keep_sample_mean=True, proj=True):
- """Estimate noise covariance matrix from raw file and events.
+def compute_covariance(epochs, keep_sample_mean=True):
+ """Estimate noise covariance matrix from epochs
The noise covariance is typically estimated on pre-stim periods
when the stim onset if defined from events.
Parameters
----------
- raw : Raw instance
- The raw data
- events : array
- The events a.k.a. the triggers
- event_ids : array-like of int
- The valid events to consider
- tmin : float
- Initial time in (s) around trigger. Ex: if tmin=0.2
- then the beginning is 200ms before trigger onset.
- tmax : float
- Final time in (s) around trigger. Ex: if tmax=0.5
- then the end is 500ms after trigger onset.
- baseline: tuple of length 2
- The time interval to apply baseline correction.
- If None do not apply it. If baseline is (a, b)
- the interval is between "a (s)" and "b (s)".
- If a is None the beginning of the data is used
- and if b is None then b is set to the end of the interval.
- If baseline is equal ot (None, None) all the time
- interval is used.
- Default is (None, None).
- reject : dict
- Rejection parameters based on peak to peak amplitude.
- Valid keys are 'grad' | 'mag' | 'eeg' | 'eog' | 'ecg'.
- If reject is None then no rejection is done.
- Values are float. Example:
- reject = dict(grad=4000e-13, # T / m (gradiometers)
- mag=4e-12, # T (magnetometers)
- eeg=40e-6, # uV (EEG channels)
- eog=250e-6 # uV (EOG channels)
- )
- flat : dict
- Rejection parameters based on flatness of signal
- Valid keys are 'grad' | 'mag' | 'eeg' | 'eog' | 'ecg'
- If flat is None then no rejection is done.
+ epochs : instance of Epochs
+ The epochs
keep_sample_mean : bool
If False data are centered at each instant before computing
the covariance.
@@ -519,33 +447,29 @@ def compute_covariance(raw, events, event_ids, tmin, tmax,
cov : instance of Covariance
The computed covariance.
"""
- # Pick all channels
- picks = pick_types(raw.info, meg=True, eeg=True, eog=True)
- if isinstance(event_ids, int):
- event_ids = list(event_ids)
- events = events.copy()
- events_numbers = events[:, 2]
- for i in event_ids[1:]:
- events_numbers[events_numbers == i] = event_ids[0]
-
data = 0.0
+ data_mean = 0.0
n_samples = 0
- epochs = Epochs(raw, events, event_ids[0], tmin, tmax, picks=picks,
- baseline=baseline, proj=proj, reject=reject, flat=flat)
+ n_epochs = 0
picks_meeg = pick_types(epochs.info, meg=True, eeg=True, eog=False)
ch_names = [epochs.ch_names[k] for k in picks_meeg]
for e in epochs:
e = e[picks_meeg]
if not keep_sample_mean:
- e -= np.mean(e, axis=0)
+ data_mean += np.sum(e, axis=0)
data += np.dot(e, e.T)
n_samples += e.shape[1]
+ n_epochs += 1
if n_samples == 0:
raise ValueError('Not enough samples to compute the noise covariance'
' matrix : %d samples' % n_samples)
- data /= n_samples - 1
+ if keep_sample_mean:
+ data /= n_samples
+ else:
+ data /= n_samples - 1
+ data -= n_samples / (1.0 - n_samples) * np.dot(data_mean, data_mean.T)
cov = Covariance(None)
cov.data = data
cov.ch_names = ch_names
diff --git a/mne/event.py b/mne/event.py
index 358c139..fefd59a 100644
--- a/mne/event.py
+++ b/mne/event.py
@@ -106,3 +106,27 @@ def find_events(raw, stim_channel='STI 014'):
idx += raw.first_samp + 1
events = np.c_[idx, np.zeros_like(idx), events_id]
return events
+
+
+def merge_events(events, ids, new_id):
+ """Merge a set of events
+
+ Parameters
+ ----------
+ events : array
+ Events
+ ids : array of int
+ The ids of events to merge
+ new_id : int
+ The new id
+
+ Returns
+ -------
+ new_events: array
+ The new events
+ """
+ events = events.copy()
+ events_numbers = events[:, 2]
+ for i in ids:
+ events_numbers[events_numbers == i] = new_id
+ return events
\ No newline at end of file
diff --git a/mne/fiff/tests/data/process_raw.sh b/mne/fiff/tests/data/process_raw.sh
index bea50e8..d7e3431 100755
--- a/mne/fiff/tests/data/process_raw.sh
+++ b/mne/fiff/tests/data/process_raw.sh
@@ -15,6 +15,10 @@ mne_process_raw --raw test_raw.fif --lowpass 40 --projoff \
mne_process_raw --raw test_raw.fif --filteroff --projon \
--savecovtag -cov --cov test.cov
+# Compute the noise covariance matrix with keepsamplemean
+mne_process_raw --raw test_raw.fif --filteroff --projon \
+ --savecovtag -km-cov --cov test_keepmean.cov
+
# Compute projection
mne_process_raw --raw test_raw.fif --events test-eve.fif --makeproj \
--projtmin -0.2 --projtmax 0.3 --saveprojtag _proj \
diff --git a/mne/fiff/tests/data/test-cov.fif b/mne/fiff/tests/data/test-cov.fif
index 7da0ed9..3806071 100755
Binary files a/mne/fiff/tests/data/test-cov.fif and b/mne/fiff/tests/data/test-cov.fif differ
diff --git a/mne/fiff/tests/data/test-km-cov.fif b/mne/fiff/tests/data/test-km-cov.fif
new file mode 100644
index 0000000..3cd3d2a
Binary files /dev/null and b/mne/fiff/tests/data/test-km-cov.fif differ
diff --git a/mne/fiff/tests/data/test.cov b/mne/fiff/tests/data/test.cov
index a9d7782..9ba2b75 100755
--- a/mne/fiff/tests/data/test.cov
+++ b/mne/fiff/tests/data/test.cov
@@ -22,7 +22,7 @@ cov {
event 1
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
def {
@@ -31,7 +31,7 @@ cov {
ignore 0
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
def {
@@ -40,7 +40,7 @@ cov {
ignore 0
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
def {
@@ -49,7 +49,7 @@ cov {
ignore 0
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
}
diff --git a/mne/fiff/tests/data/test.cov b/mne/fiff/tests/data/test_keepmean.cov
similarity index 88%
copy from mne/fiff/tests/data/test.cov
copy to mne/fiff/tests/data/test_keepmean.cov
index a9d7782..c94cda8 100755
--- a/mne/fiff/tests/data/test.cov
+++ b/mne/fiff/tests/data/test_keepmean.cov
@@ -14,6 +14,7 @@ cov {
magReject 4e-12
eegReject 80e-6
eogReject 150e-6
+ keepsamplemean
#
# What to include in the covariance matrix?
#
@@ -22,7 +23,7 @@ cov {
event 1
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
def {
@@ -31,7 +32,7 @@ cov {
ignore 0
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
def {
@@ -40,7 +41,7 @@ cov {
ignore 0
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
def {
@@ -49,7 +50,7 @@ cov {
ignore 0
tmin -0.2
tmax 0.0
- basemin -0.2
+ basemin -0.1
basemax 0
}
}
diff --git a/mne/tests/test_cov.py b/mne/tests/test_cov.py
index 8dc5337..e3638b8 100644
--- a/mne/tests/test_cov.py
+++ b/mne/tests/test_cov.py
@@ -3,12 +3,16 @@ import os.path as op
from numpy.testing import assert_array_almost_equal
from scipy import linalg
-import mne
+from .. import Covariance, read_cov, Epochs, merge_events, \
+ find_events, write_cov_file, compute_raw_data_covariance, \
+ compute_covariance
from ..fiff import fiff_open, read_evoked, Raw
from ..datasets import sample
cov_fname = op.join(op.dirname(__file__), '..', 'fiff', 'tests', 'data',
'test-cov.fif')
+cov_km_fname = op.join(op.dirname(__file__), '..', 'fiff', 'tests', 'data',
+ 'test-km-cov.fif')
raw_fname = op.join(op.dirname(__file__), '..', 'fiff', 'tests', 'data',
'test_raw.fif')
erm_cov_fname = op.join('mne', 'fiff', 'tests', 'data',
@@ -20,13 +24,13 @@ def test_io_cov():
"""
fid, tree, _ = fiff_open(cov_fname)
cov_type = 1
- cov = mne.read_cov(fid, tree, cov_type)
+ cov = read_cov(fid, tree, cov_type)
fid.close()
- mne.write_cov_file('cov.fif', cov)
+ write_cov_file('cov.fif', cov)
fid, tree, _ = fiff_open('cov.fif')
- cov2 = mne.read_cov(fid, tree, cov_type)
+ cov2 = read_cov(fid, tree, cov_type)
fid.close()
assert_array_almost_equal(cov['data'], cov2['data'])
@@ -36,8 +40,8 @@ def test_cov_estimation_on_raw_segment():
"""Estimate raw on continuous recordings (typically empty room)
"""
raw = Raw(raw_fname)
- cov = mne.compute_raw_data_covariance(raw)
- cov_mne = mne.Covariance(erm_cov_fname)
+ cov = compute_raw_data_covariance(raw)
+ cov_mne = Covariance(erm_cov_fname)
assert cov_mne.ch_names == cov.ch_names
print (linalg.norm(cov.data - cov_mne.data, ord='fro')
/ linalg.norm(cov.data, ord='fro'))
@@ -49,13 +53,23 @@ def test_cov_estimation_with_triggers():
"""Estimate raw with triggers
"""
raw = Raw(raw_fname)
- events = mne.find_events(raw)
+ events = find_events(raw)
event_ids = [1, 2, 3, 4]
- cov = mne.compute_covariance(raw, events, event_ids, tmin=-0.2, tmax=0,
- reject=dict(grad=10000e-13, mag=4e-12,
- eeg=80e-6, eog=150e-6),
- keep_sample_mean=True, proj=True)
- cov_mne = mne.Covariance(cov_fname)
+ reject = dict(grad=10000e-13, mag=4e-12, eeg=80e-6, eog=150e-6)
+
+ events = merge_events(events, event_ids, 1234)
+ epochs = Epochs(raw, events, 1234, tmin=-0.2, tmax=0,
+ baseline=(-0.2, -0.1), proj=True,
+ reject=reject)
+
+ cov = compute_covariance(epochs, keep_sample_mean=True)
+ cov_mne = Covariance(cov_km_fname)
+ assert cov_mne.ch_names == cov.ch_names
+ assert (linalg.norm(cov.data - cov_mne.data, ord='fro')
+ / linalg.norm(cov.data, ord='fro')) < 0.005
+
+ cov = compute_covariance(epochs, keep_sample_mean=False)
+ cov_mne = Covariance(cov_fname)
assert cov_mne.ch_names == cov.ch_names
assert (linalg.norm(cov.data - cov_mne.data, ord='fro')
/ linalg.norm(cov.data, ord='fro')) < 0.06
@@ -73,7 +87,7 @@ def test_whitening_cov():
# Reading
evoked = read_evoked(ave_fname, setno=0, baseline=(None, 0))
- cov = mne.Covariance(cov_fname)
+ cov = Covariance(cov_fname)
cov.get_whitener(evoked.info)
# XXX : test something
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/python-mne.git
More information about the debian-med-commit
mailing list