Python mne.morph_data函数代码示例

def test_morph_data():
    """Test morphing of data
    """
    subject_from = 'sample'
    subject_to = 'fsaverage'
    fname = op.join(data_path, 'MEG', 'sample', 'sample_audvis-meg')
    stc_from = read_source_estimate(fname)
    stc_from.crop(0.09, 0.1)  # for faster computation
    # After running this:
    #    stc_from.save('%s_audvis-meg-cropped' % subject_from)
    # this was run from a command line:
    #    mne_make_movie --stcin sample_audvis-meg-cropped-lh.stc
    #        --subject sample --morph fsaverage --smooth 12 --morphgrade 3
    #        --stc fsaverage_audvis-meg-cropped
    # XXX These files should eventually be moved to the sample dataset and
    # removed from mne/fiff/tests/data/
    fname = op.join(op.dirname(__file__), '..', 'fiff', 'tests', 'data',
                    'fsaverage_audvis-meg-cropped')
    stc_to = read_source_estimate(fname)
    stc_to1 = morph_data(subject_from, subject_to, stc_from,
                            grade=3, smooth=12, buffer_size=1000)
    stc_to1.save('%s_audvis-meg' % subject_to)
    stc_to2 = morph_data(subject_from, subject_to, stc_from,
                            grade=3, smooth=12, buffer_size=3)
    # indexing silliness here due to mne_make_movie's indexing oddities
    assert_array_almost_equal(stc_to.data, stc_to1.data[:, 0][:, None], 5)
    assert_array_almost_equal(stc_to1.data, stc_to2.data)
    # make sure precomputed morph matrices work
    vertices_to = grade_to_vertices(subject_to, grade=3)
    morph_mat = compute_morph_matrix(subject_from, subject_to,
                                     stc_from.vertno, vertices_to,
                                     smooth=12)
    stc_to3 = morph_data_precomputed(subject_from, subject_to,
                                     stc_from, vertices_to, morph_mat)
    assert_array_almost_equal(stc_to1.data, stc_to3.data)

    mean_from = stc_from.data.mean(axis=0)
    mean_to = stc_to1.data.mean(axis=0)
    assert_true(np.corrcoef(mean_to, mean_from).min() > 0.999)

    # test two types of morphing:
    # 1) make sure we can fill by morphing
    stc_to5 = morph_data(subject_from, subject_to, stc_from,
                            grade=None, smooth=12, buffer_size=3)
    assert_true(stc_to5.data.shape[0] == 163842 + 163842)

    # 2) make sure we can specify vertices
    vertices_to = [np.arange(10242), np.arange(10242)]
    stc_to3 = morph_data(subject_from, subject_to, stc_from,
                            grade=vertices_to, smooth=12, buffer_size=3)
    stc_to4 = morph_data(subject_from, subject_to, stc_from,
                            grade=5, smooth=12, buffer_size=3)
    assert_array_almost_equal(stc_to3.data, stc_to4.data)

def test_morph_data():
    """Test morphing of data
    """
    subject_from = 'sample'
    subject_to = 'fsaverage'
    fname = op.join(data_path, 'MEG', 'sample', 'sample_audvis-meg')
    stc_from = read_source_estimate(fname, subject='sample')
    fname = op.join(data_path, 'MEG', 'sample', 'fsaverage_audvis-meg')
    stc_to = read_source_estimate(fname)
    # make sure we can specify grade
    stc_from.crop(0.09, 0.1)  # for faster computation
    stc_to.crop(0.09, 0.1)  # for faster computation
    stc_to1 = stc_from.morph(subject_to, grade=3, smooth=12, buffer_size=1000,
                             subjects_dir=subjects_dir)
    stc_to1.save(op.join(tempdir, '%s_audvis-meg' % subject_to))
    # make sure we can specify vertices
    vertices_to = grade_to_vertices(subject_to, grade=3)
    stc_to2 = morph_data(subject_from, subject_to, stc_from,
                         grade=vertices_to, smooth=12, buffer_size=1000,
                         subjects_dir=subjects_dir)
    # make sure we can use different buffer_size
    stc_to3 = morph_data(subject_from, subject_to, stc_from,
                         grade=vertices_to, smooth=12, buffer_size=3,
                         subjects_dir=subjects_dir)

    assert_array_almost_equal(stc_to.data, stc_to1.data, 5)
    assert_array_almost_equal(stc_to1.data, stc_to2.data)
    assert_array_almost_equal(stc_to1.data, stc_to3.data)
    # make sure precomputed morph matrices work
    morph_mat = compute_morph_matrix(subject_from, subject_to,
                                     stc_from.vertno, vertices_to,
                                     smooth=12, subjects_dir=subjects_dir)
    stc_to3 = stc_from.morph_precomputed(subject_to, vertices_to, morph_mat)
    assert_array_almost_equal(stc_to1.data, stc_to3.data)

    mean_from = stc_from.data.mean(axis=0)
    mean_to = stc_to1.data.mean(axis=0)
    assert_true(np.corrcoef(mean_to, mean_from).min() > 0.999)

    # make sure we can fill by morphing
    stc_to5 = morph_data(subject_from, subject_to, stc_from, grade=None,
                         smooth=12, buffer_size=3, subjects_dir=subjects_dir)
    assert_true(stc_to5.data.shape[0] == 163842 + 163842)

    # test morphing to the same subject
    stc_to6 = stc_from.morph(subject_from, grade=stc_from.vertno, smooth=1,
                             subjects_dir=subjects_dir)
    mask = np.ones(stc_from.data.shape[0], dtype=np.bool)
    # XXX: there is a bug somewhere that causes a difference at 2 vertices..
    mask[6799] = False
    mask[6800] = False
    assert_array_almost_equal(stc_from.data[mask], stc_to6.data[mask], 5)

def apply_STC_epo(fnepo, event, method='MNE', snr=1.0, min_subject='fsaverage',
                  subjects_dir=None):

    from mne import morph_data
    from mne.minimum_norm import read_inverse_operator, apply_inverse_epochs

    fnlist = get_files_from_list(fnepo)
    # loop across all filenames
    for fname in fnlist:
        fn_path = os.path.split(fname)[0]
        name = os.path.basename(fname)
        subject = name.split('_')[0]
        min_dir = subjects_dir + '/%s' %min_subject
        snr = snr
        lambda2 = 1.0 / snr ** 2
        stcs_path = min_dir + '/stcs/%s/%s/' % (subject,event)
        reset_directory(stcs_path)
        # fn_inv = fname[:fname.rfind('-ave.fif')] + ',ave-inv.fif'
        fn_inv = fn_path + '/%s_epo-inv.fif' %subject

        # noise_cov = mne.read_cov(fn_cov)
        epo = mne.read_epochs(fname)
        epo.pick_types(meg=True, ref_meg=False)
        inv = read_inverse_operator(fn_inv)
        stcs = apply_inverse_epochs(epo, inv, lambda2, method,
                            pick_ori='normal')
        s = 0
        while s < len(stcs):
            stc_morph = morph_data(subject, min_subject, stcs[s])
            stc_morph.save(stcs_path + '/trial%s_fsaverage'
                           % (str(s)), ftype='stc')
            s = s + 1

def apply_inverse(fnevo, method='dSPM', snr=3.0, event='LLst', 
                  baseline=False, btmin=-0.3, btmax=-0.1, min_subject='fsaverage'):
    '''  
        Parameter
        ---------
        fnevo: string or list
            The evoked file with ECG, EOG and environmental noise free.
        method: inverse method, 'MNE' or 'dSPM'
        event: string
            The event name related with epochs.
        min_subject: string
            The subject name as the common brain.
        snr: signal to noise ratio for inverse solution. 
    '''
    #Get the default subjects_dir
    from mne.minimum_norm import apply_inverse
    fnlist = get_files_from_list(fnevo)
    # loop across all filenames
    for fname in fnlist:
        fn_path = os.path.split(fname)[0]
        name = os.path.basename(fname)
        stc_name = name[:name.rfind('-ave.fif')] 
        subject = name.split('_')[0]
        subject_path = subjects_dir + '/%s' %subject
        min_dir = subjects_dir + '/%s' %min_subject
        fn_trans = fn_path + '/%s-trans.fif' % subject
        fn_cov = fn_path + '/%s_empty,nr-cov.fif' % subject
        fn_src = subject_path + '/bem/%s-ico-5-src.fif' % subject
        fn_bem = subject_path + '/bem/%s-5120-5120-5120-bem-sol.fif' % subject
        snr = snr
        lambda2 = 1.0 / snr ** 2 
        #noise_cov = mne.read_cov(fn_cov)
        [evoked] = mne.read_evokeds(fname)
        noise_cov = mne.read_cov(fn_cov)
        # this path used for ROI definition
        stc_path = min_dir + '/%s_ROIs/%s' %(method,subject)
        #fn_cov = meg_path + '/%s_empty,fibp1-45,nr-cov.fif' % subject
        set_directory(stc_path)
        noise_cov = mne.cov.regularize(noise_cov, evoked.info,
                                        mag=0.05, grad=0.05, proj=True)
        fwd_ev = mne.make_forward_solution(evoked.info, trans=fn_trans,
                                            src=fn_src, bem=fn_bem,
                                            fname=None, meg=True, eeg=False,
                                            mindist=5.0, n_jobs=2,
                                            overwrite=True)
        fwd_ev = mne.convert_forward_solution(fwd_ev, surf_ori=True)
        forward_meg_ev = mne.pick_types_forward(fwd_ev, meg=True, eeg=False)
        inverse_operator_ev = mne.minimum_norm.make_inverse_operator(
            evoked.info, forward_meg_ev, noise_cov,
            loose=0.2, depth=0.8)
        # Compute inverse solution
        stc = apply_inverse(evoked, inverse_operator_ev, lambda2, method,
                            pick_ori=None)
        # Morph STC
        stc_morph = mne.morph_data(subject, min_subject, stc, grade=5, smooth=5)
        stc_morph.save(stc_path + '/%s' % (stc_name), ftype='stc')
        if baseline == True:
            stc_base = stc_morph.crop(btmin, btmax)
            stc_base.save(stc_path + '/%s_%s_baseline' % (subject, event), ftype='stc')

def test_morph_data():
    """Test morphing of data
    """
    subject_from = 'sample'
    subject_to = 'fsaverage'
    fname = op.join(data_path, 'MEG', 'sample', 'sample_audvis-meg')
    stc_from = SourceEstimate(fname)
    stc_from.crop(0.09, 0.1)  # for faster computation
    stc_to = morph_data(subject_from, subject_to, stc_from,
                            grade=3, smooth=12, buffer_size=1000)
    stc_to.save('%s_audvis-meg' % subject_to)

    stc_to2 = morph_data(subject_from, subject_to, stc_from,
                            grade=3, smooth=12, buffer_size=3)
    assert_array_almost_equal(stc_to.data, stc_to2.data)

    mean_from = stc_from.data.mean(axis=0)
    mean_to = stc_to.data.mean(axis=0)
    assert_true(np.corrcoef(mean_to, mean_from).min() > 0.999)

def test_morph_data():
    """Test morphing of data
    """
    subject_from = 'sample'
    subject_to = 'fsaverage'
    fname = op.join(data_path, 'MEG', 'sample', 'sample_audvis-meg')
    stc_from = read_source_estimate(fname, subject='sample')
    fname = op.join(data_path, 'MEG', 'sample', 'fsaverage_audvis-meg')
    stc_to = read_source_estimate(fname)
    # make sure we can specify grade
    stc_from.crop(0.09, 0.1)  # for faster computation
    stc_to.crop(0.09, 0.1)  # for faster computation
    stc_to1 = stc_from.morph(subject_to, grade=3, smooth=12, buffer_size=1000,
                             subjects_dir=subjects_dir)
    stc_to1.save(op.join(tempdir, '%s_audvis-meg' % subject_to))
    # make sure we can specify vertices
    vertices_to = grade_to_vertices(subject_to, grade=3)
    stc_to2 = morph_data(subject_from, subject_to, stc_from,
                         grade=vertices_to, smooth=12, buffer_size=1000,
                         subjects_dir=subjects_dir)
    # make sure we can use different buffer_size
    stc_to3 = morph_data(subject_from, subject_to, stc_from,
                         grade=vertices_to, smooth=12, buffer_size=3,
                         subjects_dir=subjects_dir)
    # indexing silliness here due to mne_make_movie's indexing oddities
    assert_array_almost_equal(stc_to.data, stc_to1.data, 5)
    assert_array_almost_equal(stc_to1.data, stc_to2.data)
    assert_array_almost_equal(stc_to1.data, stc_to3.data)
    # make sure precomputed morph matrices work
    morph_mat = compute_morph_matrix(subject_from, subject_to,
                                     stc_from.vertno, vertices_to,
                                     smooth=12, subjects_dir=subjects_dir)
    stc_to3 = stc_from.morph_precomputed(subject_to, vertices_to, morph_mat)
    assert_array_almost_equal(stc_to1.data, stc_to3.data)

    mean_from = stc_from.data.mean(axis=0)
    mean_to = stc_to1.data.mean(axis=0)
    assert_true(np.corrcoef(mean_to, mean_from).min() > 0.999)

    # make sure we can fill by morphing
    stc_to5 = morph_data(subject_from, subject_to, stc_from, grade=None,
                         smooth=12, buffer_size=3, subjects_dir=subjects_dir)
    assert_true(stc_to5.data.shape[0] == 163842 + 163842)

def _morphed_epochs_files(params):
    subject, cond_name, stc_file_name, inverse_method, subjects_dir = params
    print('morphing {}'.format(stc_file_name))
    epoch_id = utils.namebase(stc_file_name).split('_')[2]
    morphed_stc_file_name = op.join(LOCAL_ROOT_DIR, 'stc_epochs_morphed',  '{}_{}_{}_{}'.format(subject, cond_name, epoch_id, inverse_method))
    if not op.isfile('{}-stc.h5'.format(morphed_stc_file_name)):
        stc = mne.read_source_estimate(stc_file_name)
        stc_morphed = mne.morph_data(subject, 'fsaverage', stc, grade=5, smooth=20,
            subjects_dir=subjects_dir)
        stc_morphed.save(morphed_stc_file_name, ftype='h5')
    else:
        print('{} {} {} already morphed'.format(subject, cond_name, epoch_id))

def test_morph_data():
    """Test morphing of data
    """
    import mne
    subject_from = 'sample'
    subject_to = 'morph'
    fname = op.join(data_path, 'MEG', 'sample', 'sample_audvis-meg')
    stc_from = mne.SourceEstimate(fname)
    stc_to = mne.morph_data(subject_from, subject_to, stc_from, 3)

    stc_to.save('%s_audvis-meg' % subject_to)

    mean_from = stc_from.data.mean(axis=0)
    mean_to = stc_to.data.mean(axis=0)
    assert np.corrcoef(mean_to, mean_from).min() > 0.99

def smooth_ttest_results(tmin, tstep, subjects_dir, inverse_method='dSPM', n_jobs=1):
    for cond_id, cond_name in enumerate(events_id.keys()):
        for patient in get_patients():
            results_file_name = op.join(LOCAL_ROOT_DIR, 'permutation_ttest_results', '{}_{}_{}_clusters.npy'.format(patient, cond_name, inverse_method))
            if op.isfile(results_file_name):
                data = np.load(results_file_name)
                print('smoothing {} {}'.format(patient, cond_name))
                fsave_vertices = [np.arange(10242), np.arange(10242)]
                stc = _make_stc(data, fsave_vertices, tmin=tmin, tstep=tstep, subject='fsaverage')
                vertices_to = mne.grade_to_vertices('fsaverage', grade=None, subjects_dir=subjects_dir)
                print(stc.data.shape, vertices_to[0].shape)
                stc_smooth = mne.morph_data('fsaverage', 'fsaverage', stc, n_jobs=n_jobs, grade=vertices_to, subjects_dir=subjects_dir)
                stc_smooth.save(op.join(LOCAL_ROOT_DIR, 'results_for_blender', '{}_{}_{}'.format(patient, cond_name, inverse_method)), ftype='h5')
            else:
                print('no results for {} {}'.format(patient, cond_name))

def morph_STC(fn_stc, grade, template='fsaverage', event='LLst', 
              baseline=True, btmin=-0.3, btmax=0.):
    from mne import read_source_estimate, morph_data
    fnlist = get_files_from_list(fn_stc) 
    for fname in fnlist:  
        name = os.path.basename(fname)
        subject = name.split('_')[0]
        stc_name = name[:name.rfind('-ave.fif')] 
        min_dir = subjects_dir + '/%s' %template
        # this path used for ROI definition
        stc_path = min_dir + '/dSPM_ROIs/%s' %(subject)
        #fn_cov = meg_path + '/%s_empty,fibp1-45,nr-cov.fif' % subject
        set_directory(stc_path) 
        # Morph STC
        stc = read_source_estimate(fname)
        stc_morph = morph_data(subject, template, stc, grade=grade)
        stc_morph.save(stc_path + '/%s' % (stc_name), ftype='stc')
        if baseline == True:
            stc_base = stc_morph.crop(btmin, btmax)
            stc_base.save(stc_path + '/%s_%s_baseline' % (subject, event[:2]), ftype='stc') 

def morph_STC(fn_stc, grade, subjects_dir, template='fsaverage', event='LLst',
              baseline=True, btmin=-0.3, btmax=0.):
    '''
        Morph individual STC into the common brain space.

        Parameter
        ------------------------------------
        fn_stc: string or list
            The path of the individual STC.
        subjects_dir: The total bath of all the subjects.
        template: string
            The subject name as the common brain.
        event: string
            The name of event
        baseline: bool
            If true, prestimulus segment from 'btmin' to 'btmax' will be saved,
            If false, no baseline segment is saved.
        btmin, btmax: float
            If 'baseline' is True, baseline is croped using this period.

    '''
    from mne import read_source_estimate, morph_data
    fnlist = get_files_from_list(fn_stc)
    for fname in fnlist:
        name = os.path.basename(fname)
        subject = name.split('_')[0]
        stc_name = name[:name.rfind('-lh.stc')]
        min_dir = subjects_dir + '/%s' % template
        # this path used for ROI definition
        stc_path = min_dir + '/dSPM_ROIs/%s' % (subject)
        # fn_cov = meg_path + '/%s_empty,fibp1-45,nr-cov.fif' % subject
        set_directory(stc_path)
        # Morph STC
        stc = read_source_estimate(fname)
        stc_morph = morph_data(subject, template, stc, grade=grade)
        stc_morph.save(stc_path + '/%s' % (stc_name), ftype='stc')
        if baseline is True:
            stc_base = stc_morph.crop(btmin, btmax)
            stc_base.save(stc_path + '/%s_%s_baseline' % (subject, event[:2]),
                          ftype='stc')

    def __call__(self, data):

        output = deepcopy(data)

        from_surf_file = dirname(dirname(self.from_surf.surf_file))
        SUBJECTS_DIR, from_surf_name = split(from_surf_file)

        if 'lh' in self.from_surf.surf_file:  # TODO: not good, we need to check
            vertices = [arange(data.data[0].shape[0]), arange(0)]
        else:
            vertices = [arange(0), arange(data.data[0].shape[0])]

        stc = SourceEstimate(atleast_2d(data.data[0]).T, vertices=vertices,
                             tstep=0, tmin=0)
        m = morph_data(from_surf_name, self.to_surf, stc,
                       subjects_dir=SUBJECTS_DIR, grade=None,
                       smooth=self.smooth, verbose=False)

        output.data[0] = squeeze(m.data, axis=1)
        output.axis['surf'][0] = arange(m.data.shape[0])

        return output

def morph_stcs_to_fsaverage(events_id, stc_per_epoch=False, inverse_method='dSPM', subjects_dir='', n_jobs=1):
    if subjects_dir is '':
        subjects_dir = os.environ['SUBJECTS_DIR']
    for subject in get_subjects():
        for cond_name in events_id.keys():
            print('morphing {}, {}'.format(subject, cond_name))
            if not stc_per_epoch:
                morphed_stc_file_name = op.join(LOCAL_ROOT_DIR, 'stc_morphed', '{}_{}_morphed_{}'.format(subject, cond_name, inverse_method))
                if op.isfile('{}-stc.h5'.format(morphed_stc_file_name)):
                    print('{} {} already morphed'.format(subject, cond_name))
                else:
                    local_stc_file_name = op.join(LOCAL_ROOT_DIR, 'stc', '{}_{}_{}'.format(subject, cond_name, inverse_method))
                    if op.isfile('{}-stc.h5'.format(local_stc_file_name)):
                        stc = mne.read_source_estimate(local_stc_file_name)
                        stc_morphed = mne.morph_data(subject, 'fsaverage', stc, grade=5, smooth=20,
                            subjects_dir=subjects_dir)
                        stc_morphed.save(morphed_stc_file_name, ftype='h5')
                    else:
                        print("can't find stc file for {}, {}".format(subject, cond_name))
            else:
                stcs = glob.glob(op.join(LOCAL_ROOT_DIR, 'stc_epochs', '{}_{}_*_{}-stc.h5'.format(subject, cond_name, inverse_method)))
                params = [(subject, cond_name, stc_file_name, inverse_method, subjects_dir) for stc_file_name in stcs]
                utils.parallel_run(pool, _morphed_epochs_files, params, n_jobs)

def apply_inverse(fn_epo, event='LLst',ctmin=0.05, ctmax=0.25, nctmin=-0.2, nctmax=0,
                  fmin=4, fmax=8, min_subject='fsaverage', STCs=False):
    """
    Inverse evokes into source space using DICS method.
    ----------
    fn_epo : epochs of raw data.
    event_id: event id related with epochs.
    ctmin: the min time for computing CSD
    ctmax: the max time for computing CSD
    fmin: min value of the interest frequency band
    fmax: max value of the interest frequency band 
    min_subject: the subject for the common brain space.
    STCs: bool, make STCs of epochs.
    """
    from mne import Epochs, pick_types
    from mne.io import Raw
    from mne.event import make_fixed_length_events
    fnlist = get_files_from_list(fn_epo)
    # loop across all filenames
    for fname in fnlist:
        subjects_dir = os.environ['SUBJECTS_DIR']
        # extract the subject infromation from the file name
        meg_path = os.path.split(fname)[0]
        name = os.path.basename(fname)
        stc_name = name[:name.rfind('-epo.fif')] 
        subject = name.split('_')[0]
        subject_path = subjects_dir + '/%s' %subject
        min_dir = subjects_dir + '/%s' %min_subject
        fn_trans = meg_path + '/%s-trans.fif' % subject
        fn_src = subject_path + '/bem/%s-ico-4-src.fif' % subject
        fn_bem = subject_path + '/bem/%s-5120-5120-5120-bem-sol.fif' % subject
        # Make sure the target path is exist
        stc_path = min_dir + '/DICS_ROIs/%s' % subject
        set_directory(stc_path)
        # Read the MNI source space
        epochs = mne.read_epochs(fname)
        evoked = epochs.average()
        forward = mne.make_forward_solution(epochs.info, trans=fn_trans,
                                            src=fn_src, bem=fn_bem,
                                            fname=None, meg=True, eeg=False,
                                            mindist=5.0, n_jobs=2,
                                            overwrite=True)
        forward = mne.convert_forward_solution(forward, surf_ori=True)
        from mne.time_frequency import compute_epochs_csd
        from mne.beamformer import dics
        data_csd = compute_epochs_csd(epochs, mode='multitaper', tmin=ctmin, tmax=ctmax, 
                                      fmin=fmin, fmax=fmax)

        noise_csd = compute_epochs_csd(epochs, mode='multitaper', tmin=nctmin, tmax=nctmax,
                                           fmin=fmin, fmax=fmax)
                
        stc = dics(evoked, forward, noise_csd, data_csd)
        from mne import morph_data
        stc_morph = morph_data(subject, min_subject, stc, grade=4, smooth=4)
        stc_morph.save(stc_path + '/%s_%d_%d' % (stc_name, fmin, fmax), ftype='stc')
        if STCs == True:
            stcs_path = stc_path + '/STCs-%s/' %event
            reset_directory(stcs_path)
            stcs = dics(epochs, forward, noise_csd, data_csd)
            s = 0
            while s < len(stcs):
                stc_morph = mne.morph_data(subject, min_subject, stcs[s], grade=4, smooth=4)
                stc_morph.save(stcs_path + '/trial_%s'
                                % (str(s)), ftype='stc')
                s = s + 1

def morph_stc(subject_from, subject_to, stc_from_file):
    stc_from = mne.read_source_estimate(stc_from_file)
    vertices_to = [np.arange(10242), np.arange(10242)]
    stc_to = mne.morph_data(subject_from, subject_to, stc_from, n_jobs=4, grade=vertices_to)
    stc_to.save("{}_{}.stc".format(stc_from_file[:-4], subject_to))