[med-svn] [python-mne] 36/376: examples with new sample data + some docstrings
Yaroslav Halchenko
debian at onerussian.com
Fri Nov 27 17:22:02 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 9c3e90b61bc84de19b9f9c125444ac4640f703c4
Author: Alexandre Gramfort <alexandre.gramfort at inria.fr>
Date: Tue Jan 18 18:34:17 2011 -0500
examples with new sample data + some docstrings
---
examples/compute_mne_inverse.py | 9 ++-
examples/read_bem_surfaces.py | 22 ++++--
examples/read_forward.py | 3 +-
examples/read_inverse.py | 26 +++----
examples/read_stc.py | 2 +-
mne/inverse.py | 151 ++++++++++++++++++++++------------------
mne/source_space.py | 2 +-
7 files changed, 119 insertions(+), 96 deletions(-)
diff --git a/examples/compute_mne_inverse.py b/examples/compute_mne_inverse.py
index 94efae0..ed97060 100644
--- a/examples/compute_mne_inverse.py
+++ b/examples/compute_mne_inverse.py
@@ -5,20 +5,19 @@ Compute MNE inverse solution
"""
print __doc__
-from mne import fiff
+import mne
-fname_inv = 'MNE-sample-data/MEG/sample/sample_audvis-ave-7-meg-inv.fif'
+fname_inv = 'MNE-sample-data/MEG/sample/sample_audvis-meg-oct-6-meg-inv.fif'
fname_data = 'MNE-sample-data/MEG/sample/sample_audvis-ave.fif'
-# inv = fiff.read_inverse_operator(fname)
setno = 0
lambda2 = 10
dSPM = True
-res = fiff.compute_inverse(fname_data, setno, fname_inv, lambda2, dSPM)
+res = mne.compute_inverse(fname_data, setno, fname_inv, lambda2, dSPM)
import pylab as pl
pl.plot(res['sol'][::100,:].T)
-pl.xlabel('time (s)')
+pl.xlabel('time (ms)')
pl.ylabel('Source amplitude')
pl.show()
diff --git a/examples/read_bem_surfaces.py b/examples/read_bem_surfaces.py
index 127f9c4..c4faa5d 100644
--- a/examples/read_bem_surfaces.py
+++ b/examples/read_bem_surfaces.py
@@ -5,19 +5,27 @@ Reading BEM surfaces from a forward solution
"""
print __doc__
-from mne import fiff
+import mne
-fname = 'MNE-sample-data/subjects/sample/bem/sample-5120-bem-sol.fif'
+fname = 'MNE-sample-data/subjects/sample/bem/sample-5120-5120-5120-bem-sol.fif'
-surf = fiff.read_bem_surfaces(fname)
+surfaces = mne.read_bem_surfaces(fname)
-print "Number of surfaces : %d" % len(surf)
+print "Number of surfaces : %d" % len(surfaces)
###############################################################################
# Show result
+head_col = (0.95, 0.83, 0.83) # light pink
+skull_col = (0.91, 0.89, 0.67)
+brain_col = (0.67, 0.89, 0.91) # light blue
+colors = [head_col, skull_col, brain_col]
+
# 3D source space
-points = surf[0]['rr']
-faces = surf[0]['tris']
from enthought.mayavi import mlab
-mlab.triangular_mesh(points[:,0], points[:,1], points[:,2], faces)
+mlab.clf()
+for c, surf in zip(colors, surfaces):
+ points = surf['rr']
+ faces = surf['tris']
+ mlab.triangular_mesh(points[:,0], points[:,1], points[:,2], faces,
+ color=c, opacity=0.3)
diff --git a/examples/read_forward.py b/examples/read_forward.py
index 75a9bad..c6377f0 100644
--- a/examples/read_forward.py
+++ b/examples/read_forward.py
@@ -7,8 +7,7 @@ print __doc__
import mne
-# fname = 'MNE-sample-data/MEG/sample/sample_audvis-ave-7-fwd.fif'
-fname = 'sm01a5-ave-oct-6-fwd.fif'
+fname = 'MNE-sample-data/MEG/sample/sample_audvis-meg-oct-6-fwd.fif'
data = mne.read_forward_solution(fname)
leadfield = data['sol']['data']
diff --git a/examples/read_inverse.py b/examples/read_inverse.py
index 96e8c47..00025ed 100644
--- a/examples/read_inverse.py
+++ b/examples/read_inverse.py
@@ -7,7 +7,7 @@ print __doc__
import mne
-fname = 'MNE-sample-data/MEG/sample/sample_audvis-ave-7-meg-inv.fif'
+fname = 'MNE-sample-data/MEG/sample/sample_audvis-meg-oct-6-meg-inv.fif'
inv = mne.read_inverse_operator(fname)
@@ -16,14 +16,16 @@ print "fMRI prior: %s" % inv['fmri_prior']
print "Number of sources: %s" % inv['nsource']
print "Number of channels: %s" % inv['nchan']
-# ###############################################################################
-# # Show result
-#
-# # 3D source space
-# lh_points = inv['src'][0]['rr']
-# lh_faces = inv['src'][0]['use_tris']
-# rh_points = inv['src'][1]['rr']
-# rh_faces = inv['src'][1]['use_tris']
-# from enthought.mayavi import mlab
-# mlab.triangular_mesh(lh_points[:,0], lh_points[:,1], lh_points[:,2], lh_faces)
-# mlab.triangular_mesh(rh_points[:,0], rh_points[:,1], rh_points[:,2], rh_faces)
+###############################################################################
+# Show result
+
+# 3D source space
+import numpy as np
+lh_points = inv['src'][0]['rr']
+lh_faces = inv['src'][0]['use_tris']
+rh_points = inv['src'][1]['rr']
+rh_faces = inv['src'][1]['use_tris']
+from enthought.mayavi import mlab
+mlab.triangular_mesh(lh_points[:,0], lh_points[:,1], lh_points[:,2], lh_faces)
+mlab.triangular_mesh(rh_points[:,0], rh_points[:,1], rh_points[:,2], rh_faces)
+mlab.show()
diff --git a/examples/read_stc.py b/examples/read_stc.py
index 17cd649..d5acd0e 100644
--- a/examples/read_stc.py
+++ b/examples/read_stc.py
@@ -10,7 +10,7 @@ print __doc__
import mne
-fname = 'MNE-sample-data/MEG/sample/sample_audvis-ave-7-meg-lh.stc'
+fname = 'MNE-sample-data/MEG/sample/sample_audvis-meg-lh.stc'
stc = mne.read_stc(fname)
diff --git a/mne/inverse.py b/mne/inverse.py
index 8e67bae..53d92b7 100644
--- a/mne/inverse.py
+++ b/mne/inverse.py
@@ -66,14 +66,14 @@ def read_inverse_operator(fname):
raise ValueError, 'Modalities not found'
inv = dict()
- inv['methods'] = tag.data;
+ inv['methods'] = tag.data
tag = find_tag(fid, invs, FIFF.FIFF_MNE_SOURCE_ORIENTATION)
if tag is None:
fid.close()
raise ValueError, 'Source orientation constraints not found'
- inv['source_ori'] = tag.data;
+ inv['source_ori'] = tag.data
tag = find_tag(fid, invs, FIFF.FIFF_MNE_SOURCE_SPACE_NPOINTS)
if tag is None:
@@ -110,26 +110,29 @@ def read_inverse_operator(fname):
fid.close()
raise ValueError, 'Singular values not found'
- inv['sing'] = tag.data
+ inv['sing'] = tag.data
inv['nchan'] = len(inv['sing'])
#
# The eigenleads and eigenfields
#
inv['eigen_leads_weighted'] = False
try:
- inv['eigen_leads'] = _read_named_matrix(fid, invs, FIFF.FIFF_MNE_INVERSE_LEADS)
+ inv['eigen_leads'] = _read_named_matrix(fid, invs,
+ FIFF.FIFF_MNE_INVERSE_LEADS)
except:
- inv['eigen_leads_weighted'] = True
- try:
- inv.eigen_leads = _read_named_matrix(fid,invs,FIFF.FIFF_MNE_INVERSE_LEADS_WEIGHTED);
- except Exception as inst:
- raise ValueError, '%s' % inst
+ inv['eigen_leads_weighted'] = True
+ try:
+ inv.eigen_leads = _read_named_matrix(fid, invs,
+ FIFF.FIFF_MNE_INVERSE_LEADS_WEIGHTED)
+ except Exception as inst:
+ raise ValueError, '%s' % inst
#
# Having the eigenleads as columns is better for the inverse calculations
#
inv['eigen_leads'] = _transpose_named_matrix(inv['eigen_leads'])
try:
- inv['eigen_fields'] = _read_named_matrix(fid, invs, FIFF.FIFF_MNE_INVERSE_FIELDS)
+ inv['eigen_fields'] = _read_named_matrix(fid, invs,
+ FIFF.FIFF_MNE_INVERSE_FIELDS)
except Exception as inst:
raise ValueError, '%s' % inst
@@ -164,13 +167,13 @@ def read_inverse_operator(fname):
FIFF.FIFFV_MNE_DEPTH_PRIOR_COV)
print '\tDepth priors read.'
except:
- inv['depth_prior'] = [];
+ inv['depth_prior'] = []
try:
inv['fmri_prior'] = read_cov(fid, invs, FIFF.FIFFV_MNE_FMRI_PRIOR_COV)
print '\tfMRI priors read.'
except:
- inv['fmri_prior'] = [];
+ inv['fmri_prior'] = []
#
# Read the source spaces
@@ -182,7 +185,7 @@ def read_inverse_operator(fname):
raise ValueError, 'Could not read the source spaces (%s)' % inst
for s in inv['src']:
- s['id'] = find_source_space_hemi(s)
+ s['id'] = find_source_space_hemi(s)
#
# Get the MRI <-> head coordinate transformation
@@ -192,14 +195,15 @@ def read_inverse_operator(fname):
fid.close()
raise ValueError, 'MRI/head coordinate transformation not found'
else:
- mri_head_t = tag.data;
+ mri_head_t = tag.data
if mri_head_t['from_'] != FIFF.FIFFV_COORD_MRI or \
mri_head_t['to'] != FIFF.FIFFV_COORD_HEAD:
mri_head_t = _invert_transform(mri_head_t)
if mri_head_t['from_'] != FIFF.FIFFV_COORD_MRI or \
mri_head_t['to'] != FIFF.FIFFV_COORD_HEAD:
fid.close()
- raise ValueError, 'MRI/head coordinate transformation not found'
+ raise ValueError, ('MRI/head coordinate transformation '
+ 'not found')
inv['mri_head_t'] = mri_head_t
#
@@ -215,7 +219,7 @@ def read_inverse_operator(fname):
#
# Number of averages is initially one
#
- inv['nave'] = 1;
+ inv['nave'] = 1
#
# We also need the SSP operator
#
@@ -223,24 +227,25 @@ def read_inverse_operator(fname):
#
# Some empty fields to be filled in later
#
- inv['proj'] = [] # This is the projector to apply to the data
- inv['whitener'] = [] # This whitens the data
- inv['reginv'] = [] # This the diagonal matrix implementing
+ inv['proj'] = [] # This is the projector to apply to the data
+ inv['whitener'] = [] # This whitens the data
+ inv['reginv'] = [] # This the diagonal matrix implementing
# regularization and the inverse
inv['noisenorm'] = [] # These are the noise-normalization factors
#
nuse = 0
for k in range(len(inv['src'])):
- try:
- inv['src'][k] = _transform_source_space_to(inv['src'][k],
+ try:
+ inv['src'][k] = _transform_source_space_to(inv['src'][k],
inv['coord_frame'], mri_head_t)
- except Exception as inst:
- fid.close()
- raise ValueError, 'Could not transform source space (%s)', inst
+ except Exception as inst:
+ fid.close()
+ raise ValueError, 'Could not transform source space (%s)', inst
- nuse += inv['src'][k]['nuse']
+ nuse += inv['src'][k]['nuse']
- print '\tSource spaces transformed to the inverse solution coordinate frame'
+ print ('\tSource spaces transformed to the inverse solution '
+ 'coordinate frame')
#
# Done!
#
@@ -267,7 +272,7 @@ def combine_xyz(vec):
raise ValueError, ('Input must be a 1D vector with '
'3N entries')
- s = _block_diag(vec[None,:], 3)
+ s = _block_diag(vec[None, :], 3)
comb = (s * s.T).diagonal()
return comb
@@ -295,12 +300,12 @@ def prepare_inverse_operator(orig, nave, lambda2, dSPM):
# Scale some of the stuff
#
scale = float(inv['nave']) / nave
- inv['noise_cov']['data'] = scale * inv['noise_cov']['data']
- inv['noise_cov']['eig'] = scale * inv['noise_cov']['eig']
+ inv['noise_cov']['data'] = scale * inv['noise_cov']['data']
+ inv['noise_cov']['eig'] = scale * inv['noise_cov']['eig']
inv['source_cov']['data'] = scale * inv['source_cov']['data']
#
if inv['eigen_leads_weighted']:
- inv['eigen_leads']['data'] = sqrt(scale) * inv['eigen_leads']['data']
+ inv['eigen_leads']['data'] = sqrt(scale) * inv['eigen_leads']['data']
print ('\tScaled noise and source covariance from nave = %d to '
@@ -309,7 +314,7 @@ def prepare_inverse_operator(orig, nave, lambda2, dSPM):
#
# Create the diagonal matrix for computing the regularized inverse
#
- inv['reginv'] = inv['sing'] / (inv['sing'] * inv['sing'] + lambda2);
+ inv['reginv'] = inv['sing'] / (inv['sing'] * inv['sing'] + lambda2)
print '\tCreated the regularized inverter'
#
# Create the projection operator
@@ -330,22 +335,23 @@ def prepare_inverse_operator(orig, nave, lambda2, dSPM):
#
nnzero = 0
for k in range(ncomp, inv['noise_cov']['dim']):
- if inv['noise_cov']['eig'][k] > 0:
- inv['whitener'][k,k] = 1.0 / sqrt(inv['noise_cov']['eig'][k])
- nnzero += 1
+ if inv['noise_cov']['eig'][k] > 0:
+ inv['whitener'][k, k] = 1.0 / sqrt(inv['noise_cov']['eig'][k])
+ nnzero += 1
#
# Rows of eigvec are the eigenvectors
#
inv['whitener'] = np.dot(inv['whitener'], inv['noise_cov']['eigvec'])
print ('\tCreated the whitener using a full noise covariance matrix '
- '(%d small eigenvalues omitted)' % (inv['noise_cov']['dim'] - nnzero))
+ '(%d small eigenvalues omitted)' % (inv['noise_cov']['dim']
+ - nnzero))
else:
#
# No need to omit the zeroes due to projection
#
for k in range(inv['noise_cov']['dim']):
- inv['whitener'][k,k] = 1.0 / sqrt(inv['noise_cov']['data'][k])
+ inv['whitener'][k, k] = 1.0 / sqrt(inv['noise_cov']['data'][k])
print ('\tCreated the whitener using a diagonal noise covariance '
'matrix (%d small eigenvalues discarded)' % ncomp)
@@ -358,12 +364,13 @@ def prepare_inverse_operator(orig, nave, lambda2, dSPM):
noise_norm = np.zeros(inv['eigen_leads']['nrow'])
if inv['eigen_leads_weighted']:
for k in range(inv['eigen_leads']['nrow']):
- one = inv['eigen_leads']['data'][k,:] * inv['reginv']
+ one = inv['eigen_leads']['data'][k, :] * inv['reginv']
noise_norm[k] = np.sum(one**2)
else:
for k in range(inv['eigen_leads']['nrow']):
one = sqrt(inv['source_cov']['data'][k]) * \
- np.sum(inv['eigen_leads']['data'][k,:] * inv['reginv'])
+ np.sum(inv['eigen_leads']['data'][k, :]
+ * inv['reginv'])
noise_norm[k] = np.sum(one**2)
#
@@ -388,31 +395,39 @@ def prepare_inverse_operator(orig, nave, lambda2, dSPM):
inv['noisenorm'] = np.diag(1.0 / np.abs(noise_norm)) # XXX
print '[done]'
else:
- inv['noisenorm'] = [];
+ inv['noisenorm'] = []
return inv
def compute_inverse(fname_data, setno, fname_inv, lambda2, dSPM=True,
nave=None):
- """
- %
- % [res] = mne_ex_compute_inverse(fname_data,setno,fname_inv,nave,lambda2,dSPM)
- %
- % An example on how to compute a L2-norm inverse solution
- % Actual code using these principles might be different because
- % the inverse operator is often reused across data sets.
- %
- %
- % fname_data - Name of the data file
- % setno - Data set number
- % fname_inv - Inverse operator file name
- % nave - Number of averages (scales the noise covariance)
- % If negative, the number of averages in the data will be
- % used
- % lambda2 - The regularization factor
- % dSPM - do dSPM?
- %
+ """Compute inverse solution
+
+ Computes a L2-norm inverse solution
+ Actual code using these principles might be different because
+ the inverse operator is often reused across data sets.
+
+ Parameters
+ ----------
+ fname: string
+ File name of the data file
+ setno: int
+ Data set number
+ fname_inv: string
+ File name of the inverse operator
+ nave: int
+ Number of averages (scales the noise covariance)
+ If negative, the number of averages in the data will be used XXX
+ lambda2: float
+ The regularization parameter
+ dSPM: bool
+ do dSPM ?
+
+ Returns
+ -------
+ res: dict
+ Inverse solution
"""
#
@@ -453,17 +468,17 @@ def compute_inverse(fname_data, setno, fname_inv, lambda2, dSPM=True,
# with the weights computed above
#
if inv['eigen_leads_weighted']:
- #
- # R^0.5 has been already factored in
- #
- print '(eigenleads already weighted)...',
- sol = np.dot(inv['eigen_leads']['data'], trans)
+ #
+ # R^0.5 has been already factored in
+ #
+ print '(eigenleads already weighted)...',
+ sol = np.dot(inv['eigen_leads']['data'], trans)
else:
- #
- # R^0.5 has to factored in
- #
- print '(eigenleads need to be weighted)...',
- sol = np.sqrt(inv['source_cov']['data'])[:,None] * \
+ #
+ # R^0.5 has to factored in
+ #
+ print '(eigenleads need to be weighted)...',
+ sol = np.sqrt(inv['source_cov']['data'])[:, None] * \
np.dot(inv['eigen_leads']['data'], trans)
@@ -471,7 +486,7 @@ def compute_inverse(fname_data, setno, fname_inv, lambda2, dSPM=True,
print 'combining the current components...',
sol1 = np.zeros((sol.shape[0]/3, sol.shape[1]))
for k in range(sol.shape[1]):
- sol1[:,k] = np.sqrt(combine_xyz(sol[:,k]))
+ sol1[:, k] = np.sqrt(combine_xyz(sol[:, k]))
sol = sol1
diff --git a/mne/source_space.py b/mne/source_space.py
index 1b10c8d..074c82a 100644
--- a/mne/source_space.py
+++ b/mne/source_space.py
@@ -144,7 +144,7 @@ def _read_one_source_space(fid, this, open_here):
fid.close()
raise ValueError, 'Vertex data not found'
- res['rr'] = tag.data
+ res['rr'] = tag.data.astype(np.float) # make it double precision for mayavi
if res['rr'].shape[0] != res['np']:
if open_here:
fid.close()
--
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