Wednesday, May 28, 2014

Day 7: Combining source spaces

Notebook

Today we'll attempt to combine the subcortical volume source space with a cortical surface source space.
The code below is last week's code with a few modifications. Namely, I'm adding the cortical source space
and then computing the forward model.

In [2]:
import numpy as np
import mne
from mne.datasets import spm_face
from mne.minimum_norm import make_inverse_operator, apply_inverse_epochs
import matplotlib.pyplot as plt
%matplotlib inline

# supress text output
mne.set_log_level(False)

# set the data directories
data_path = spm_face.data_path()
subjects_dir = data_path + '/subjects'
subject = 'spm'
subject_dir = data_path + '/subjects/' + subject
meg_dir = data_path + '/MEG/spm'

# get the data files
bem_fname = subject_dir + '/bem/spm-5120-5120-5120-bem-sol.fif'
mri_fname = meg_dir + '/SPM_CTF_MEG_example_faces1_3D_raw-trans.fif'
epo_fname = meg_dir + '/SPM_CTF_MEG_example_faces1_3D_epo.fif'

# load the epoch data
epochs = mne.read_epochs(epo_fname)

# get positions and orientations of the right amygdala (subcortical)
pos = mne.source_space.get_segment_positions('spm', 'Right-Amygdala',
                                             random_ori=True,
                                             subjects_dir=subjects_dir)

# setup the right amygdala volume source space
vol_src = mne.setup_volume_source_space('spm', pos=pos)

# setup the cortical surface source space
src = mne.setup_source_space('spm', overwrite=True)

# combine the source spaces
src.append(vol_src[0])

# setup the forward model
forward = mne.make_forward_solution(epochs.info, mri=mri_fname, src=src,
                                    bem=bem_fname)
forward = mne.convert_forward_solution(forward, surf_ori=True)

# Compute inverse solution
snr = 5.0
lambda2 = 1.0 / snr ** 2
method = 'dSPM'

# estimate noise covarariance
noise_cov = mne.compute_covariance(epochs.crop(None, 0, copy=True))

# make the inverse operator
inverse_operator = make_inverse_operator(epochs.info, forward, noise_cov,
                                         loose=0.2, depth=0.8)

# Apply inverse solution to both stim types
stc_faces = apply_inverse_epochs(epochs['faces'], inverse_operator, lambda2,
                                 method)
stc_scrambled = apply_inverse_epochs(epochs['scrambled'], inverse_operator,
                                     lambda2, method)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-2-03275dc2df3d> in <module>()
     57 # Apply inverse solution to both stim types
     58 stc_faces = apply_inverse_epochs(epochs['faces'], inverse_operator, lambda2,
---> 59                                  method)
     60 stc_scrambled = apply_inverse_epochs(epochs['scrambled'], inverse_operator,
     61                                      lambda2, method)

/Users/Alan/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/mne/minimum_norm/inverse.pyc in apply_inverse_epochs(epochs, inverse_operator, lambda2, method, label, nave, pick_ori, return_generator, verbose, pick_normal)

/Users/Alan/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/mne/utils.pyc in verbose(function, *args, **kwargs)
    390         return ret
    391     else:
--> 392         ret = function(*args, **kwargs)
    393         return ret
    394 

/Users/Alan/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/mne/minimum_norm/inverse.pyc in apply_inverse_epochs(epochs, inverse_operator, lambda2, method, label, nave, pick_ori, return_generator, verbose, pick_normal)
   1000     if not return_generator:
   1001         # return a list
-> 1002         stcs = [stc for stc in stcs]
   1003 
   1004     return stcs

/Users/Alan/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/mne/minimum_norm/inverse.pyc in _apply_inverse_epochs_gen(epochs, inverse_operator, lambda2, method, label, nave, pick_ori, verbose, pick_normal)
    945 
    946         stc = _make_stc(sol, vertices=vertno, tmin=tmin, tstep=tstep,
--> 947                         subject=subject)
    948 
    949         yield stc

/Users/Alan/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/mne/source_estimate.pyc in _make_stc(data, vertices, tmin, tstep, subject)
    344                                 tstep=tstep, subject=subject)
    345     else:
--> 346         raise ValueError('vertices has to be either a list with one or two '
    347                          'arrays or an array')
    348     return stc

ValueError: vertices has to be either a list with one or two arrays or an array

As you can see, the program crashed.

Digging deeper into the mne code, it looks like the function is setup to read either a volume or a surface
source space, but not both. Fortunately, another mne contributor shared his processing stream for
combining source spaces. In his implementation, he combines the forward operators, not the source spaces.

More on that implementation tomorrow.

No comments:

Post a Comment