Plot external source analysis

💬 Support & Discussions
,
1
  • MNE-Python version: 0.23.1
  • operating system: Ubuntu 18.04.5 LTS as server running Jupyter and macbook BigSur 11.5 as client (with firefox 91.0.2).

Dear MNE community,

After performing source localization in my data, I extracted the time series per label using stc.in_label() and performed some analysis outside MNE, ending up with 68 values (one per region), stored in the some array external_computations. Now I would like to plot these values in source space. What would be the easiest way to do that in MNE? I couldn’t find anything similar in the tutorials or in the discussions here.

I did hack something up after inspecting the mne.SourceEstimate class, which seems to be working, but it just looks… wrong :slight_smile:

Here is the snippet that seems to do what I want, but I hope there is a better way of doing that:

# Load original source estimates
labels = mne.read_labels_from_annot(subject, parc='aparc', hemi='both', subjects_dir=freesurfer_path)
stclh = read_source_estimate('./sources_sLORETA-lh.stc')
stcrh = read_source_estimate('./sources_sLORETA-rh.stc')

# Crop to get only one sample
stclh = stclh.crop(0, .0001)
stcrh = stcrh.crop(0, .0001)

# Set original source estimate values to zero
stclh.data[:] = 0
stcrh.data[:] = 0

# Update the values in each dipole of each region with the external computations in "external_computations", which is an array with 68 values, one per region
 
for (il, l) in enumerate(labels):
    if (l.name[-2:] == 'lh'):
        lvert = stclh.in_label(l).vertices[0]
        inter = np.intersect1d(lvert, stclh.vertices[0], return_indices=True)
        stclh.data[inter[2], 0] = external_computations[il]
        
    if (l.name[-2:] == 'rh'):
        lvert = stcrh.in_label(l).vertices[1]
        inter = np.intersect1d(lvert, stcrh.vertices[1], return_indices=True)
        stcrh.data[inter[2]+len(stclh.vertices[0]), 0] = external_computations[il]

Now in theory I can do the plots with stclh.plot() and stcrh.plot(). But I’m not sure this is correct, and it looks overly complicated. For instance I have to sum the number of vertices in left hemisphere (len(stclh.vertices[0])) to the positions identified in the right hemisphere (stcrh.in_label(l).vertices[1]) before indexing the actual data (stcrh.data[]). Is this the correct way to access the data of a particular region in the right hemisphere?

Best regards,
Leonardo

2

I don’t think there is a much simpler solution.

Alex

1 Like
3

Ok, I was just not sure if this was right, for instance I got the fact that the vertices in the right hemisphere were offset simply by trial and error (just indexing without the shift gave me “christmas tree” activity). I guess I’ll stick to this solution then. Tks!