MNE: Close-up on forward solutions

In previous notebooks we computed the forward solution operator given the raw data, transformed data, BEM model and source space. In this notebook we investigate deeper the forward solution operator. Let us import the forward solution operator we previously computed.

import mne

subjects_dir = r'C:\Users\hz3752\PycharmProjects\mne_bids_pipeline\data\meg\Sub-0037\sub-0037-fwd.fif'

fwd = mne.read_forward_solution(subjects_dir, include=(), exclude=(), ordered=True, verbose=None)

Reading forward solution from C:\Users\hz3752\PycharmProjects\mne_bids_pipeline\data\meg\Sub-0037\sub-0037-fwd.fif...
    Reading a source space...
    Computing patch statistics...
    Patch information added...
    Distance information added...
    [done]
    Reading a source space...
    Computing patch statistics...
    Patch information added...
    Distance information added...
    [done]
    2 source spaces read
    Desired named matrix (kind = 3523) not available
    Read MEG forward solution (5124 sources, 207 channels, free orientations)
    Source spaces transformed to the forward solution coordinate frame

Two possible source orientations exist in MNE, documented in source_ori in fwd:

• FIFF.FIFFV_MNE_FIXED_ORI: Fixed orientation mode, the orientation of the source currents is constrained to be perpendicular to the cortical surface. This means that the direction of the current is assumed to be known a priori and is fixed during the inverse solution process. The primary advantage of this approach is its simplicity and reduced computational load. By constraining the direction of the sources, the number of variables in the inverse problem is effectively halved, which can make the solutions more stable under conditions of noisy data or when the data are limited.

• When to use fixed orientation? Fixed orientation is often used when a strong prior knowledge of the orientation of the cortical currents exists (such as orientation aligned with the local sulcal or gyral geometry) or when computational simplicity is required.

• FIFF.FIFFV_MNE_FREE_ORI: Free Orientation: In contrast, free orientation allows the orientation of the source currents to vary and be estimated as part of the source localization process. This mode does not impose any constraints on the direction of the currents, thereby enabling the estimation of the most probable orientation based on the data itself. This can potentially lead to more accurate reconstructions of the underlying neural activity, particularly in complex scenarios where the orientation of the sources cannot be easily predicted.

• When to use free orientation?: Free orientation is particularly useful in research settings where the precise modeling of neural activity is crucial, and where there is enough high-quality data to support the increased complexity of the inverse solution. For the data we red, let us print which source orientation was used:

print(fwd['source_ori'])

2 (FIFFV_MNE_FREE_ORI)

The location of the sources are in source_rr

print(fwd['source_rr'])

[[-0.0210264  -0.06770435  0.06570835]
 [-0.01708252 -0.06216683  0.07929588]
 [-0.01315619 -0.06514898  0.07275735]
 ...
 [ 0.00288066  0.02885441  0.04047859]
 [ 0.01054338  0.0624362   0.0432792 ]
 [ 0.01025996  0.08508232  0.04505657]]

The forward solution operator are in sol they are of shape (n_channels, n_sources * n_orientation) where n_orientation is the number of possible orientations of the source (3 orientations for the free-orientation setup and 1 for the fixed orientation setup), let us print the shape of the solution

print(fwd['sol']['data'].shape)

(207, 15372)

This means that we used 15372/3 = 5124 sources and the solutions are computed for 207 MEG gradiometers

import mne
# Load your raw data
raw = mne.io.read_raw_fif(r'C:\Users\hz3752\PycharmProjects\mne_bids_pipeline\data\meg\Sub-0037\sub-01_01-eyes-closed-raw.fif', preload=True)

Opening raw data file C:\Users\hz3752\PycharmProjects\mne_bids_pipeline\data\meg\Sub-0037\sub-01_01-eyes-closed-raw.fif...
    Range : 0 ... 321999 =      0.000 ...   321.999 secs
Ready.
Reading 0 ... 321999  =      0.000 ...   321.999 secs...

#Add covariance matrix computation

#cov = mne.compute_covariance

cov = mne.compute_raw_covariance(raw)
cov.save(r'C:\Users\hz3752\PycharmProjects\mne_bids_pipeline\data\meg\Sub-0037\sub-01-cov.fif')

Using up to 1610 segments
Number of samples used : 322000
[done]

 inv = mne.minimum_norm.make_inverse_operator(raw.info, fwd, cov, depth=None, loose='auto', fixed=False) #fixed=False: Ignoring dipole direction.

Converting forward solution to surface orientation
    Average patch normals will be employed in the rotation to the local surface coordinates....
    Converting to surface-based source orientations...
    [done]
Computing inverse operator with 207 channels.
    207 out of 207 channels remain after picking
Selected 207 channels
Applying loose dipole orientations to surface source spaces: 0.2
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 8.2e-13 (2.2e-16 eps * 207 dim * 18  max singular value)
    Estimated rank (mag): 207
    MAG: rank 207 computed from 207 data channels with 0 projectors
    Setting small MAG eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 14.1305
    scaling factor to adjust the trace = 7.29474e+27 (nchan = 207 nzero = 3)

subj='sub-01'
        #--------------------------STCs--------------------------------#

print '%s: Creating STCs...'%subj
        os.makedirs('STC/%s' %subj)
        for ev in evoked:
            stc = mne.minimum_norm.apply_inverse(ev, inv, lambda2=lambda2, method='dSPM')
            # mophing stcs to the fsaverage using precomputed matrix method:
            vertices_to = mne.grade_to_vertices('fsaverage', grade=4, 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.morph_data_precomputed(subject_from=subj, subject_to='fsaverage', stc_from=stc, vertices_to=vertices_to, morph_mat=morph_mat)
            stc_morph.save('STC/%s/%s_%s_dSPM' %(subj,subj,ev.comment))
            del stc, stc_morph
        print '>> DONE CREATING STCS FOR SUBJ=%s'%subj
        print '-----------------------------------------\n'

        #deleting variables
        del epochs_rej, evoked, info, trans, src, fwd, cov, inv