Hi all,
Is there any formula to find out the surface area corresponding to a dipole within a source space? When I am looking at the MNE documentation regarding ico-5 source space, its written as:
Source space: āico5ā
Sources per hemisphere: 10242
Source spacing: 3.1 mm
Surface area per source : 9.8 mm2 (is there any general formula to calculate this?)
Thanks,
There is this info:
>>> import mne
>>> src = mne.read_source_spaces(mne.datasets.sample.data_path() / 'subjects' / 'fsaverage' / 'bem' / 'fsaverage-ico-5-src.fif', patch_stats=True)
>>> src[0].keys()
dict_keys(['id', 'type', 'np', 'ntri', 'coord_frame', 'rr', 'nn', 'tris', 'nuse', 'inuse', 'vertno', 'nuse_tri', 'use_tris', 'nearest', 'nearest_dist', 'pinfo', 'patch_inds', 'dist', 'dist_limit', 'subject_his_id', 'tri_area', 'tri_cent', 'tri_nn', 'use_tri_cent', 'use_tri_nn', 'use_tri_area'])
>>> src[0]["use_tri_area"].shape
(20480,)
>>> src[0]["use_tris"].shape
(20480, 3)
So this could get you the surface area per triangle in the decimated space. But you want it per vertex not per triangle.
Really the āsurface area per sourceā can be thought of as ātotal surface area / total number of sourcesā. So looking at it that way, you could calculate it as:
sum(s["tri_area"].sum() for s in src) / sum(s['nuse'] for s in src) * 1e6
6.367790682029399
The * 1e6 is a conversion from āsources per square mā to sources per square mm". For an oct-6 for example youād get:
>>> src = mne.read_source_spaces(mne.datasets.sample.data_path() / 'subjects' / 'sample' / 'bem' / 'sample-oct-6-src.fif', patch_stats=True)
>>> sum(s["use_tri_area"].sum() for s in src) / sum(s['nuse'] for s in src) * 1e6
20.69581503362706
Iām not sure why the value for the ico-5 doesnāt match, assuming that manual entry is for fsaverage (because it will vary by subject). Even that for sample wouldnāt be 9.8:
>>> sum(s["use_tri_area"].sum() for s in src) / 20484 * 1e6
8.280750830677961
ā¦ oh in those last two I errantly used use_tri_area ā if we use tri_area we do get a value that matches (if you floor instead of round):
>>> sum(s["tri_area"].sum() for s in src) / 20484 * 1e6
9.882712190552494
@larsoner Thanks a lot for the explanation.
It works with sample data but not the 'fsaverge'
In [38]: src= read_source_spaces('/home/dip_meg/mne_data/MNE-sample-data/subjects/fsaverage/bem/fsaverage-ico-5-src.fif')
Reading a source space...
[done]
Reading a source space...
[done]
2 source spaces read
In [39]: src
Out[39]: <SourceSpaces: [<surface (lh), n_vertices=163842, n_used=10242>, <surface (rh), n_vertices=163842, n_used=10242>] MRI (surface RAS) coords, subject 'fsaverage', ~21.9 MB>
In [40]: sum(s["tri_area"].sum() for s in src) / 20484 * 1e6
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
Cell In[40], line 1
----> 1 sum(s["tri_area"].sum() for s in src) / 20484 * 1e6
Cell In[40], line 1, in <genexpr>(.0)
----> 1 sum(s["tri_area"].sum() for s in src) / 20484 * 1e6
KeyError: 'tri_area'
In [41]: sum(s["use_tri_area"].sum() for s in src) / 20484 * 1e6
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
Cell In[41], line 1
----> 1 sum(s["use_tri_area"].sum() for s in src) / 20484 * 1e6
Cell In[41], line 1, in <genexpr>(.0)
----> 1 sum(s["use_tri_area"].sum() for s in src) / 20484 * 1e6
KeyError: 'use_tri_area'
It seems like use_tri_area or tri_area attributes are missing in the src. Am I missing something?
thanks!
You have to pass patch_stats=True to the read_source_spaces call
I can print the values now but for ico-5 for fsaverage, it doesnāt give the same value of 9.8 as documented in MNE. Would that mean the value is rather based on the sample subject src and not for fsaverage?
In [3]: from mne import read_source_spaces
In [4]: src= read_source_spaces('/home/dip_meg/mne_data/MNE-fsaverage-data/fsave
...: rage/bem/fsaverage-ico-5-src.fif', patch_stats=True)
Reading a source space...
[done]
Completing triangulation info...
[done]
Completing selection triangulation info...
[done]
Reading a source space...
[done]
Completing triangulation info...
[done]
Completing selection triangulation info...
[done]
2 source spaces read
In [5]: src
Out[5]: <SourceSpaces: [<surface (lh), n_vertices=163842, n_used=10242>, <surface (rh), n_vertices=163842, n_used=10242>] MRI (surface RAS) coords, subject 'fsaverage', ~59.1 MB>
In [6]: src.info
Out[6]:
{'fname': '/home/dip_meg/mne_data/MNE-fsaverage-data/fsaverage/bem/fsaverage-ico-5-src.fif',
'working_dir': '/autofs/space/megmix_002/users/mluessi/data/code/mne-scripts',
'command_line': 'mne_make_source_space --surf lh.white:rh.white --ico 5 --src /cluster/fusion/mluessi/MNE-sample-data-build/MNE-sample-data/subjects/fsaverage/bem/fsaverage-ico-5-src.fif'}
In [7]: sum(s["tri_area"].sum() for s in src) / 20484 * 1e6
Out[7]: 6.367790682029399
In [8]: sum(s["use_tri_area"].sum() for s in src) / 20484 * 1e6
Out[8]: 6.246565416014217
I guess it must be. Itāll depend on the brain surface area since the number of vertices for ico-5 is always the sameā¦
Yes, you are right. I looked it up with my subjects data.
thanks for helping.
@larsoner sorry, I have a follow up question.
If I look at the mne.setup_source_space ā MNE 1.5.0 documentation , MNE creates (at-least with default setup) source space at the interface between gray and white matter, i.e. at the white surface (decimated surface, more suitable for morphological feature)? and not at the middle surface, i.e. a surface that runs at the mid-distance between white and pial??
I am looking at the mne code to compute Cortical surface area for the triangulation.
If I understand the code properly, then this means (I put the comments at the side):
def _complete_source_space_info(this, verbose=None):
"""Add more info on surface."""
# Main triangulation
logger.info(" Completing triangulation info...")
this["tri_area"] = np.zeros(this["ntri"])
r1 = this["rr"][this["tris"][:, 0], :]
r2 = this["rr"][this["tris"][:, 1], :]
r3 = this["rr"][this["tris"][:, 2], :]
this["tri_cent"] = (r1 + r2 + r3) / 3.0
this["tri_nn"] = fast_cross_3d((r2 - r1), (r3 - r1))
this["tri_area"] = _normalize_vectors(this["tri_nn"]) / 2.0
logger.info("[done]")
# Selected triangles
logger.info(" Completing selection triangulation info...")
if this["nuse_tri"] > 0:
r1 = this["rr"][this["use_tris"][:, 0], :] >>>>>>>>> a triangular face say, r1r2r3; vertex coordinates
r2 = this["rr"][this["use_tris"][:, 1], :] >>>>>>>>> vertex coordinates
r3 = this["rr"][this["use_tris"][:, 2], :] >>>>>>>>> vertex coordinates
this["use_tri_cent"] = (r1 + r2 + r3) / 3.0 >>>>>>>>> Conversion from facewise to vertexwise each vertex one-third of the sum of the areas of all faces that meet at that vertex
this["use_tri_nn"] = fast_cross_3d((r2 - r1), (r3 - r1)) >>>>>>>>>>> computes the triangle vector normal
this["use_tri_area"] = np.linalg.norm(this["use_tri_nn"], axis=1) / 2.0 >>>>>>>>>>> computes The triangle areas from vector normal in m2
logger.info("[done]")
Am I correct? I am following this paper : https://s3.us-east-2.amazonaws.com/brainder/publications/2012/winkler2012_facewise_surface_area.pdf to understand the logic.
Highly appreciate your help.
Dip
That sounds right, the white surface should be the white/gray interface
I donāt have the bandwidth to ācheck your workā so-to-speak on this ā itās probably better to convince yourself using more reading and/or checking against other libraries (like vtk) anyway!
Thanks! I will do that.