Day 29: Downsampling the volume source

Today I'm trying to figure out how to subsample a subcortical volume. In previous posts, the number of vertices making
up the amygdala (a relatively small structure) is about the same as the number of vertices making up the entire cortex.
Ideally, I'd like to restrict the number of vertices in a subcortical volume to the same scale as the cortical vertices.

Currently, mne.setup_volume_source space can fill a surface with equally spaced points. However, it first converts the
surface to a sphere. In the following example, the amygdala is not spherical.

In [25]:

%matplotlib inline
import numpy as np
import nibabel as nib
from mne.datasets import spm_face

# get the aseg file
data_path = spm_face.data_path()
aseg_fname = data_path + '/subjects/spm/mri/aseg.mgz'

# read aseg file using nibabel
aseg = nib.load(aseg_fname)
hdr = aseg.get_header()
dat = aseg.get_data()

# get indices to region of interest
ix = dat == 18
X, Y, Z = np.where(ix)

In [14]:

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# plot the region of interest
ax = plt.axes(projection='3d')
ax.scatter(x, y, z)

Out[14]:

<mpl_toolkits.mplot3d.art3d.Patch3DCollection at 0x47e9590>

One idea I had was to reduce the number of vertices was to exclude vertices who's neighbors are also vertices.
This doesn't solve any minimum spacing issues, but its a start. Unfortunately, the number of vertices doesn't even
drop by half.

In [28]:

ex, ey, ez = [], [], []
for i, x in enumerate(X):
    y = Y[i]
    z = Z[i]
    
    s = ix[x-1:x+2, y-1:y+2, z-1:z+2].sum()
    if s < 27:
        ex.append(x)
        ey.append(y)
        ez.append(z)

In [30]:

ax = plt.axes(projection='3d')
ax.scatter(ex, ey, ez)

Out[30]:

<mpl_toolkits.mplot3d.art3d.Patch3DCollection at 0x5da00d0>

In [36]:

len(X), len(ex)

Out[36]:

(1840, 1187)

In [37]:

hdr

Out[37]:

<nibabel.freesurfer.mghformat.MGHHeader at 0x47f2b50>

In [38]:

hdr.get_vox2ras_tkr()

Out[38]:

array([[  -1.,    0.,    0.,  128.],
       [   0.,    0.,    1., -128.],
       [   0.,   -1.,    0.,  128.],
       [   0.,    0.,    0.,    1.]], dtype=float32)

My thought for next week is to create a 3D matrix of equally spaced points with the same 256x256x256 dimensions
and look for overlap with the amygdala.