[med-svn] [python-mne] 312/376: ENH : make SSP computation work + adding tests

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 398b4c835b3f8c2878bf94bd9190aae468821de9
Author: Alexandre Gramfort <alexandre.gramfort at inria.fr>
Date:   Thu Jul 7 11:47:22 2011 +0200

    ENH : make SSP computation work + adding tests
---
 mne/fiff/proj.py                   | 17 ++++++---
 mne/fiff/tests/data/process_raw.sh |  2 +-
 mne/fiff/tests/test_proj.py        | 76 ++++++++++++++++++++++++--------------
 3 files changed, 60 insertions(+), 35 deletions(-)

diff --git a/mne/fiff/proj.py b/mne/fiff/proj.py
index bf2b8bb..9cb31ca 100644
--- a/mne/fiff/proj.py
+++ b/mne/fiff/proj.py
@@ -302,8 +302,7 @@ def compute_spatial_vectors(epochs, n_grad=2, n_mag=2, n_eeg=2):
     projs: list
         List of projection vectors
     """
-    data = epochs.get_data()
-    data = data.swapaxes(0, 1).reshape(data.shape[1], -1)
+    data = sum(np.dot(e, e.T) for e in epochs)  # compute data covariance
 
     mag_ind = pick_types(epochs.info, meg='mag')
     grad_ind = pick_types(epochs.info, meg='grad')
@@ -322,18 +321,24 @@ def compute_spatial_vectors(epochs, n_grad=2, n_mag=2, n_eeg=2):
     grad_names, mag_names, eeg_names = ([epochs.ch_names[k] for k in ind]
                                      for ind in [grad_ind, mag_ind, eeg_ind])
 
+    event_id = epochs.event_id
     projs = []
-    for n, ind, names in zip([n_grad, n_mag, n_eeg],
+    for n, ind, names, desc in zip([n_grad, n_mag, n_eeg],
                       [grad_ind, mag_ind, eeg_ind],
-                      [grad_names, mag_names, eeg_names]):
+                      [grad_names, mag_names, eeg_names],
+                      ['planar', 'axial', 'eeg']):
         if n == 0:
             continue
-        U = linalg.svd(data[ind], full_matrices=False,
+        data_ind = data[ind][:,ind]
+        U = linalg.svd(data_ind, full_matrices=False,
                                          overwrite_a=True)[0][:, :n]
         for k, u in enumerate(U.T):
             proj_data = dict(col_names=names, row_names=None,
                              data=u[np.newaxis, :], nrow=1, ncol=u.size)
-            proj = dict(active=True, data=proj_data, desc='MEG %s' % k, kind=1)
+            this_desc = "%s-%-d-%-.3f-%-.3f-PCA-%02d" % (desc, event_id,
+                                             epochs.tmin, epochs.tmax, k + 1)
+            print "Adding projection: %s" % this_desc
+            proj = dict(active=True, data=proj_data, desc=this_desc, kind=1)
             projs.append(proj)
 
     return projs
diff --git a/mne/fiff/tests/data/process_raw.sh b/mne/fiff/tests/data/process_raw.sh
index d7e3431..a13a078 100755
--- a/mne/fiff/tests/data/process_raw.sh
+++ b/mne/fiff/tests/data/process_raw.sh
@@ -23,4 +23,4 @@ mne_process_raw --raw test_raw.fif --filteroff --projon \
 mne_process_raw --raw test_raw.fif --events test-eve.fif --makeproj \
            --projtmin -0.2 --projtmax 0.3 --saveprojtag _proj \
            --projnmag 1 --projngrad 1 --projevent 1 \
-           --projmagrej 6000 --projgradrej 5000
+           --projmagrej 600000 --projgradrej 500000  --filteroff
diff --git a/mne/fiff/tests/test_proj.py b/mne/fiff/tests/test_proj.py
index 972df0c..53d8deb 100644
--- a/mne/fiff/tests/test_proj.py
+++ b/mne/fiff/tests/test_proj.py
@@ -1,9 +1,10 @@
 import os.path as op
 
+import numpy as np
 from numpy.testing import assert_array_almost_equal
 
 from .. import Raw, pick_types, compute_spatial_vectors
-from ..proj import make_projector
+from ..proj import make_projector, read_proj
 from ..open import fiff_open
 from ... import read_events, Epochs
 
@@ -11,30 +12,49 @@ raw_fname = op.join(op.dirname(__file__), 'data', 'test_raw.fif')
 event_fname = op.join(op.dirname(__file__), 'data', 'test-eve.fif')
 proj_fname = op.join(op.dirname(__file__), 'data', 'test_proj.fif')
 
-# XXX
-# def test_compute_spatial_vectors():
-#     """Test SSP computation
-#     """
-#     event_id, tmin, tmax = 1, -0.2, 0.3
-# 
-#     raw = Raw(raw_fname)
-#     events = read_events(event_fname)
-#     exclude = ['MEG 2443', 'EEG 053']
-#     picks = pick_types(raw.info, meg=True, eeg=False, stim=False, eog=False,
-#                             exclude=exclude)
-#     epochs = Epochs(raw, events, event_id, tmin, tmax, picks=picks,
-#                         baseline=(None, 0), proj=False,
-#                         reject=dict(mag=5000e-15, grad=16e-10))
-# 
-#     projs = compute_spatial_vectors(epochs, n_grad=1, n_mag=2, n_eeg=2)
-# 
-#     proj, nproj, U = make_projector(projs, epochs.ch_names, bads=[])
-#     assert nproj == 3
-#     assert U.shape[1] == 3
-# 
-#     epochs.info['projs'] += projs
-#     evoked = epochs.average()
-#     evoked.save('foo.fif')
-# 
-#     fid, tree, _ = fiff_open(proj_fname)
-#     projs = read_proj(fid, tree)
+
+def test_compute_spatial_vectors():
+    """Test SSP computation"""
+    event_id, tmin, tmax = 1, -0.2, 0.3
+
+    raw = Raw(raw_fname)
+    events = read_events(event_fname)
+    exclude = []
+    bad_ch = 'MEG 2443'
+    picks = pick_types(raw.info, meg=True, eeg=False, stim=False, eog=False,
+                            exclude=exclude)
+    epochs = Epochs(raw, events, event_id, tmin, tmax, picks=picks,
+                        baseline=None, proj=False)
+
+    projs = compute_spatial_vectors(epochs, n_grad=1, n_mag=1, n_eeg=0)
+
+    fid, tree, _ = fiff_open(proj_fname)
+    projs2 = read_proj(fid, tree)
+
+    for k, (p1, p2) in enumerate(zip(projs, projs2)):
+        assert p1['desc'] == p2['desc']
+        assert p1['data']['col_names'] == p2['data']['col_names']
+        # compare with sign invariance
+        p1_data = p1['data']['data'] * np.sign(p1['data']['data'][0, 0])
+        p2_data = p2['data']['data'] * np.sign(p2['data']['data'][0, 0])
+        if bad_ch in p1['data']['col_names']:
+            bad = p1['data']['col_names'].index('MEG 2443')
+            mask = np.ones(p1_data.size, dtype=np.bool)
+            mask[bad] = False
+            p1_data = p1_data[:, mask]
+            p2_data = p2_data[:, mask]
+        corr = np.corrcoef(p1_data, p2_data)[0, 1]
+        assert_array_almost_equal(corr, 1.0, 7)
+
+    # test that you can compute the projection matrix
+    proj, nproj, U = make_projector(projs, epochs.ch_names, bads=[])
+    assert nproj == 2
+    assert U.shape[1] == 2
+
+    # test that you can save them
+    epochs.info['projs'] += projs
+    evoked = epochs.average()
+    evoked.save('foo.fif')
+
+    fid, tree, _ = fiff_open(proj_fname)
+    projs = read_proj(fid, tree)

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