Misalignment coregging BrainVision CapTrak with fsaverage

Hello all,

Iā€™m trying to align EEG sensor positions, as digitized by a BrainVision CapTrak system, with fsaverage in a coregistration stage. Iā€™m using the function mne.channels.read_dig_captrak to read in the digitized sensor positions, and saving a .fif file of the EEG data with the re-aligned montage. Then, Iā€™m then aligning the resulting .fif file with fsaverage in the MNE coreg GUI. Everything looks fine, but when I make ico-4 source space and use the trans file to align the source and sensor space to make the STCs, the electrode positions seem to be shifted along the z-axis.

Hereā€™s a screenshot of what I see in the coreg GUI, and a screenshot of what I see when I use mne.viz.plot_alignment:


Oddly, if I open the mne coreg GUI, DONā€™T align the sensors, and then save a trans file, I get much better alignment.

Hereā€™s a screenshot of what happens if I donā€™t bother trying to align the CapTrak sensors with the average brain, and the projected alignment when I use mne.viz.plot_alignment. This seems much more reasonable? (Am I introducing an unnecessary step?)

Any clue? Could it be an issue with the fiducial points, or an issue with the MRI and EEG coordinates?

Iā€™m uploading code for aligning the CapTrak coordinates with the EEG montage and the source reconstruction code.

Hereā€™s our script that relates the BrainVision CapTrak files (.bvct) with the .fif files produced by our preprocessing pipeline:

import mne
import main
import os

subjDir = os.path.join(main.ROOT, 'MRI')
mne.utils.set_config('SUBJECTS_DIR', subjDir, set_env=True)
mne.set_config('SUBJECTS_DIR', subjDir, set_env=True)

for subject in main.subjects:
	#raw = mne.io.read_raw_brainvision(os.path.join('EEG',subject,subject+'.vhdr'))
	raw = mne.read_epochs(os.path.join('EEG',subject,subject+'-ica-good-epo.fif'))
	montage = mne.channels.read_dig_captrak(os.path.join('EEG',subject,subject+'.bvct'))
	raw.set_montage(montage)
	raw.save(os.path.join('EEG',subject,subject+'-ica-good-epo_CapTrak.fif'), overwrite=True)
	layout = raw.plot_sensors(show_names=True)
	if not os.path.exists(os.path.join('Plots','Individuals',subject)):
		os.mkdir(os.path.join('Plots','Individuals',subject))
	layout.savefig(os.path.join('Plots','Individuals',subject,subject+'_layout.png'))

print('Please run <mne coreg> from the terminal, and coreg the subjects!')
print('Make sure to save <subject_name>-trans.fif and to save the MRI file')

And hereā€™s the source reconstruction script:

import mne#, #eelbrain, 
import os, glob, pickle
import numpy as np
import pandas as pd
from os.path import join

import main

mne.set_log_level(verbose='WARNING')

#=========Edit here=========#

SNR = 3 # 3 for ANOVAs, 2 for regressions
fixed = False # False for orientation free (=unsigned), True for fixed orientation (=signed)

os.chdir(main.ROOT) #setting current dir

# '''==========================================================================='''
# '''                             PART 3: Create STCs                           '''
# '''==========================================================================='''

#---------------------------average STCs------------------------------------#

subjects_dir = os.path.join(main.ROOT, 'MRI')

for subject in main.subjects:
    if os.path.exists(os.path.join('STC_GrAvg', subject)):
        print('Looks like we started subject %s already! Lets just check in :)' %subject)
    else:
        os.makedirs(os.path.join('STC_GrAvg', subject))

        print(">> STCs for subj=%s:"%subject)
        print('Importing data...')

        raw = mne.io.read_raw_brainvision(os.path.join('EEG',subject,subject+'.vhdr'))
        info = raw.info

        epochs_ica = mne.read_epochs(os.path.join('EEG', subject, subject+'-ica-good-epo_CapTrak.fif'))

        covEpochs = epochs_ica.copy()

        covEpochs = covEpochs.crop(tmin=main.epoch_baseline[0],tmax=main.epoch_baseline[1])
#        covEpochs.reject = dict(grad=2e-12, mag=2e-2)

        trans = mne.read_trans(os.path.join('EEG', subject, subject+'-trans.fif'))

        bem_fname = os.path.join(subjects_dir, subject, 'bem', subject+'-5120-5120-5120-bem-sol.fif')
        src_fname = os.path.join(subjects_dir, subject, 'bem', subject+'-ico-4-src.fif')
        fwd_fname = os.path.join('EEG', subject, subject+'-fwd.fif')
        cov_fname = os.path.join('EEG', subject, subject+'-cov.fif')

         #---------------------------get evoked------------------------------------#
        print('%s: Creating evoked responses' %subject)
        evoked = []
        conditions = main.event_id.keys()
        for cond in conditions:
            evoked.append(epochs_ica[cond].average())
        print('Done.')

        #----------------------Source space---------------------------#
        print('Generating source space...')
        if os.path.isfile(src_fname):
            print('src for subj = %s already exists, loading file...' %subject)
            src = mne.read_source_spaces(fname=src_fname)
            print('Done.')
        else:
            print('src for subj = %s does not exist, creating file...' %subject)
            src = mne.setup_source_space(subject=subject, spacing='ico4', subjects_dir=subjects_dir)
            src.save(src_fname, overwrite=True)
            print('Done. File saved.')  

        #-------------------------- BEM ------------------------------#
        if not os.path.isfile(bem_fname):
#        if True:
            print('BEM for subj=%s does not exists, creating...' %subject)
            conductivity = (0.3, 0.006, 0.3) # for three layers layer
            model = mne.make_bem_model(subject=subject, ico=4, conductivity=conductivity, subjects_dir=subjects_dir)
            bem = mne.make_bem_solution(model)
            mne.write_bem_solution(bem_fname, bem, overwrite=True)

        #--------------------Forward solution-------------------------#
        print('Creating forward solution...')
        if os.path.isfile(fwd_fname):
            print('forward solution for subj=%s exists, loading file.' %subject)
            fwd = mne.read_forward_solution(fwd_fname)
            print('Done.')
        else:
            print('forward solution for subj=%s does not exist, creating file.' %subject)
            fwd = mne.make_forward_solution(info=info, trans=trans, src=src, bem=bem_fname, ignore_ref=True)
            mne.write_forward_solution(fwd_fname, fwd)
            print('Done. File saved.')


        #----------------------Covariance------------------------------#
        print('Getting covariance')
        if os.path.isfile(cov_fname):
            print('covariance matrix for subj=%s exists, loading file...' %subject)
            cov = mne.read_cov(cov_fname)

            evokedPlot = epochs_ica.average()
            whitePlot = mne.viz.plot_evoked_white(evokedPlot, cov, show=False)
            whitePlot.savefig(os.path.join('Plots','Individuals',subject, subject+'-whitenedCov.png'))


#            evoked.plot_white(cov)

            print('Done.')
        else:
            print('covariance matrix for subj=%s does not exist, creating file...' %subject)
            cov = mne.compute_covariance(covEpochs,main.epoch_baseline[0],tmax=main.epoch_baseline[1], method=['shrunk', 'diagonal_fixed', 'empirical'])
            cov.save(cov_fname)

            evokedPlot = epochs_ica.average()
            whitePlot = mne.viz.plot_evoked_white(evokedPlot, cov, show=False)
            whitePlot.savefig(os.path.join('Plots','Individuals',subject, subject+'-whitenedCov.png'))


#            evoked.plot_white(cov)            

            
            print('Done. File saved.')


        #---------------------Inverse operator-------------------------#
        print('Getting inverse operator')
        if fixed == True:
            fwd = mne.convert_forward_solution(fwd, surf_ori=True)

        inv = mne.minimum_norm.make_inverse_operator(info, fwd, cov, depth=None, #loose=None, 
          fixed=fixed) #fixed=False: Ignoring dipole direction.
        lambda2 = 1.0 / SNR ** 2.0

        #--------------------------STCs--------------------------------#

        print('%s: Creating STCs...'%subject)

#        os.makedirs('STCAvg/%s' %subj)
        for ev in evoked:

            ev.set_eeg_reference(projection=True)

            stc = mne.minimum_norm.apply_inverse(ev, inv, lambda2=lambda2, 
            #method='sLORETA'#, 
            method='sLORETA',
#                pick_ori = "normal"
            ) # This will get you signed data
            # mophing stcs to the fsaverage using precomputed matrix method:

            # Miriam's 
#            vertices_to = mne.grade_to_vertices('fsaverage', grade = [np.arange(2562), np.arange(2562)], subjects_dir=subjects_dir) #fsaverage's source space
#            morph = mne.compute_source_morph(stc, subject_from=subj, subject_to='fsaverage', subjects_dir=subjects_dir, spacing = [np.arange(2562), np.arange(2562)],
#                src_to = vertices_to
#                )
#            stc_morph = morph.apply(stc)

            # Graham's
            vertices_to = mne.grade_to_vertices('fsaverage', grade=4, subjects_dir=subjects_dir) #fsaverage's source space
            morph_mat = mne.compute_source_morph(stc, subject_from = subject, subject_to = 'fsaverage', subjects_dir=subjects_dir, spacing = 4)
            stc_morph = morph_mat.apply(stc)


            # Julien's
#            vertices_to = mne.grade_to_vertices('fsaverage', grade=5, subjects_dir=subjects_dir) #fsaverage's source space
#            morph_mat = mne.compute_morph_matrix(subject_from=subj, subject_to='fsaverage', vertices_from=stc.vertices, vertices_to=vertices_to, subjects_dir=subjects_dir)
#            stc_morph = mne.compute_source_morph(stc, subject_from=subj, subject_to='fsaverage', #stc_from=stc, 
                #vertices_to=vertices_to, 
                #morph_mat=morph_mat
#                spacing=4,
#                subjects_dir = subjects_dir
#                ).apply(stc)

            # Dustin's
#mne.compute_source_morph(stc, subject_from=subj, )

            stc_morph.save(os.path.join('STC_GrAvg', subject,'%s_%s_dSPM' %(subject,ev.comment)))
#            del stc, stc_morph
        print('>> DONE CREATING STCS FOR SUBJ=%s'%subject)
        print('-----------------------------------------\n')

        #deleting variables
        del evoked, info, trans, src, fwd, cov, inv
  • MNE Version 1.4
  • Mac OS 12.6

Iā€™m sorry! I uploading the wrong photo. The second photo is what happens when I plot the source/sensor space alignment with correctly coregistered EEG/MRIs:

As you can see, something seems a bit off :slight_smile:

The problem is solved.

The issue is that my forward solution uses the info file produced by the raw object, i.e., the EEG with the standard montage. Instead, it should have been the info file from the object epochs_ica. Epochs_ica is the object that has the digitized electrode positions :man_facepalming:

2 Likes