How to interpret a volumetric brain segmentation¶
In a previous lesson with digital images, we explored the notion of a map/parcellation. In this chapter we’ll do the same for NIfTI data, and see how parcellations are implemented in this context. Let’s begin by loading a parcellation generated by FreeSurfer (Fischl 2012) a software package for neuroimaging data which can, among other things, automatedly segment structural brain data (i.e. a T1) into known structures and/or areas. Note that this parcellation is specific to a particular subject and data source. In this case, it was derived from the subject and source (T1) that we viewed in the previous lesson.
#this code ensures that we can navigate the WiMSE repo across multiple systems
import subprocess
import os
#get top directory path of the current git repository, under the presumption that
#the notebook was launched from within the repo directory
gitRepoPath=subprocess.check_output(['git', 'rev-parse', '--show-toplevel']).decode('ascii').strip()
#move to the top of the directory
os.chdir(gitRepoPath)
import nibabel as nib
#establish path to parcellation
atlasPath=os.path.join(gitRepoPath,'exampleData','parc.nii.gz')
#load it as an object
atlasImg = nib.load(atlasPath)
#extract the header information
atlasHeader = atlasImg.header
print(atlasHeader)
<class 'nibabel.nifti1.Nifti1Header'> object, endian='<'
sizeof_hdr : 348
data_type : b''
db_name : b''
extents : 0
session_error : 0
regular : b''
dim_info : 0
dim : [ 3 256 256 256 1 1 1 1]
intent_p1 : 0.0
intent_p2 : 0.0
intent_p3 : 0.0
intent_code : none
datatype : int32
bitpix : 32
slice_start : 0
pixdim : [-1. 1. 1. 1. 1.8 1. 1. 1. ]
vox_offset : 0.0
scl_slope : nan
scl_inter : nan
slice_end : 0
slice_code : unknown
xyzt_units : 10
cal_max : 0.0
cal_min : 0.0
slice_duration : 0.0
toffset : 0.0
glmax : 0
glmin : 0
descrip : b'FreeSurfer Jan 18 2017'
aux_file : b''
qform_code : scanner
sform_code : scanner
quatern_b : 0.0
quatern_c : 0.70710677
quatern_d : -0.70710677
qoffset_x : 127.0
qoffset_y : -145.0
qoffset_z : 147.0
srow_x : [ -1. 0. 0. 127.]
srow_y : [ 0. 0. 1. -145.]
srow_z : [ 0. -1. 0. 147.]
intent_name : b''
magic : b'n+1'
Let’s take a quick moment to look at how this compares to our source T1.
#this code ensures that we can navigate the WiMSE repo across multiple systems
import subprocess
import os
#get top directory path of the current git repository, under the presumption that
#the notebook was launched from within the repo directory
gitRepoPath=subprocess.check_output(['git', 'rev-parse', '--show-toplevel']).decode('ascii').strip()
#move to the top of the directory
os.chdir(gitRepoPath)
#set path to T1
t1Path=os.path.join(gitRepoPath,'exampleData','t1.nii.gz')
#load the t1
t1img = nib.load(t1Path)
#extract the header info
T1header = t1img.header
#print the output
print(T1header)
<class 'nibabel.nifti1.Nifti1Header'> object, endian='<'
sizeof_hdr : 348
data_type : b''
db_name : b''
extents : 0
session_error : 0
regular : b'r'
dim_info : 0
dim : [ 3 182 218 182 1 1 1 1]
intent_p1 : 0.0
intent_p2 : 0.0
intent_p3 : 0.0
intent_code : none
datatype : int16
bitpix : 16
slice_start : 0
pixdim : [-1. 1. 1. 1. 1.8 0. 0. 0. ]
vox_offset : 0.0
scl_slope : nan
scl_inter : nan
slice_end : 0
slice_code : unknown
xyzt_units : 10
cal_max : 0.0
cal_min : 0.0
slice_duration : 0.0
toffset : 0.0
glmax : 0
glmin : 0
descrip : b'FSL5.0'
aux_file : b''
qform_code : mni
sform_code : mni
quatern_b : 0.0
quatern_c : 1.0
quatern_d : 0.0
qoffset_x : 90.0
qoffset_y : -126.0
qoffset_z : -72.0
srow_x : [-1. 0. 0. 90.]
srow_y : [ 0. 1. 0. -126.]
srow_z : [ 0. 0. 1. -72.]
intent_name : b''
magic : b'n+1'
Let’s compare a few of these fields. Specifically:
Data dimensions (dim)
Voxel dimensions (pixdim)
Affine (qoffset_ & srow_)
#get parcellation data shape and print it
print('Atlas Data Dimensions')
atlasDataDimensions=atlasImg.shape
print(atlasDataDimensions)
#get T1 data shape and print it
print('T1 Data Dimensions')
t1DataDimensions=t1img.shape
print(t1DataDimensions)
print('')
#print parcellation voxel dimensions
print('Atlas Voxel dimensions (in mm)')
parcellationVoxelDims=atlasImg.header.get_zooms()
print(parcellationVoxelDims)
#print T1 voxel dimensions
T1VoxelDims=t1img.header.get_zooms()
print('T1 Voxel dimensions (in mm)')
print(T1VoxelDims)
print('')
#display volume occupied for Atlas
print('Volume of space represented by parcellation')
print(f'{atlasDataDimensions[0]*parcellationVoxelDims[0]} mm by {atlasDataDimensions[1]*parcellationVoxelDims[1]} mm by {atlasDataDimensions[2]*parcellationVoxelDims[2]} mm' )
#display volume occupied for T1
print('Volume of space represented by T1')
print(f'{t1DataDimensions[0]*T1VoxelDims[0]} mm by {t1DataDimensions[1]*T1VoxelDims[1]} mm by {t1DataDimensions[2]*T1VoxelDims[2]} mm' )
Atlas Data Dimensions
(256, 256, 256)
T1 Data Dimensions
(182, 218, 182)
Atlas Voxel dimensions (in mm)
(1.0, 1.0, 1.0)
T1 Voxel dimensions (in mm)
(1.0, 1.0, 1.0)
Volume of space represented by parcellation
256.0 mm by 256.0 mm by 256.0 mm
Volume of space represented by T1
182.0 mm by 218.0 mm by 182.0 mm
Above, we see that there are obvious discrepancies between the internal structures of these data objects. Although both are composed of isometric, 1 mm voxels, only the parcellation forms a proper cube (with all sides being of the same length/size). Furthermore, the parcellation characterizes a larger volume of space.
We’ll come back to this as we try and align these images. For now, let’s move towards a consideration of the actual data content of the parcellation data object.
import numpy as np
#open it as a numpy data array
atlasData = atlasImg.get_fdata()
#set the print option so it isn't printing in scientific notation
np.set_printoptions(suppress=True)
#condense to unique values
uniqueAtlasEntries=np.unique(atlasData).astype(int)
#print the output
print('Numbers corresponding to unique labels in this atlas')
print('')
print(uniqueAtlasEntries)
Numbers corresponding to unique labels in this atlas
[ 0 2 4 5 7 8 10 11 12 13 14 15
16 17 18 24 26 28 30 31 41 43 44 46
47 49 50 51 52 53 54 58 60 62 63 77
85 251 252 253 254 255 1000 2000 11101 11102 11103 11104
11105 11106 11107 11108 11109 11110 11111 11112 11113 11114 11115 11116
11117 11118 11119 11120 11121 11122 11123 11124 11125 11126 11127 11128
11129 11130 11131 11132 11133 11134 11135 11136 11137 11138 11139 11140
11141 11143 11144 11145 11146 11147 11148 11149 11150 11151 11152 11153
11154 11155 11156 11157 11158 11159 11160 11161 11162 11163 11164 11165
11166 11167 11168 11169 11170 11171 11172 11173 11174 11175 12101 12102
12103 12104 12105 12106 12107 12108 12109 12110 12111 12112 12113 12114
12115 12116 12117 12118 12119 12120 12121 12122 12123 12124 12125 12126
12127 12128 12129 12130 12131 12132 12133 12134 12135 12136 12137 12138
12139 12140 12141 12143 12144 12145 12146 12147 12148 12149 12150 12151
12152 12153 12154 12155 12156 12157 12158 12159 12160 12161 12162 12163
12164 12165 12166 12167 12168 12169 12170 12171 12172 12173 12174 12175]
Just as we did before with our T1 data, let’s also consider the frequency of these values as well.
#transform data array for use
unwrappedAtlasData=np.ndarray.flatten(atlasData)
#this time we'll use pandas to do our count, as this will help with plotting
import pandas as pd
#convert single vector numpy array to dataframe
AtlasLabelsDataframe=pd.DataFrame(unwrappedAtlasData)
#count unique values via pandas function
AtlasLabelCounts=AtlasLabelsDataframe.iloc[:,0].value_counts().reset_index()
#rename the columns
AtlasLabelCounts=AtlasLabelCounts.rename(columns={"index":"Label Number",0:"voxel count"})
#round, for the sake of clarity, and resort them by label number
AtlasLabelCounts=AtlasLabelCounts.astype(np.int64).sort_values(by="Label Number").reset_index(drop=True)
#use itables for interactive tables,
#NOTE SEEMS TO BE BROKEN, FIX AT LATER DATE
import itables
#itables.show(AtlasLabelCounts,columnDefs=[{"width": "10px", "targets": [0,1]}])
#itables.show(AtlasLabelCounts)
AtlasLabelCounts.head(25)
| Label Number | voxel count | |
|---|---|---|
| 0 | 0 | 15382362 |
| 1 | 2 | 268665 |
| 2 | 4 | 7431 |
| 3 | 5 | 280 |
| 4 | 7 | 16651 |
| 5 | 8 | 66521 |
| 6 | 10 | 9147 |
| 7 | 11 | 4001 |
| 8 | 12 | 5313 |
| 9 | 13 | 2597 |
| 10 | 14 | 824 |
| 11 | 15 | 1403 |
| 12 | 16 | 23399 |
| 13 | 17 | 5490 |
| 14 | 18 | 2022 |
| 15 | 24 | 1077 |
| 16 | 26 | 466 |
| 17 | 28 | 4934 |
| 18 | 30 | 12 |
| 19 | 31 | 593 |
| 20 | 41 | 265864 |
| 21 | 43 | 7156 |
| 22 | 44 | 145 |
| 23 | 46 | 16910 |
| 24 | 47 | 65814 |
| 25 | 49 | 8647 |
| 26 | 50 | 4075 |
| 27 | 51 | 5558 |
| 28 | 52 | 2353 |
| 29 | 53 | 5268 |
| 30 | 54 | 2205 |
| 31 | 58 | 573 |
| 32 | 60 | 4873 |
| 33 | 62 | 5 |
| 34 | 63 | 596 |
| 35 | 77 | 556 |
| 36 | 85 | 164 |
| 37 | 251 | 1147 |
| 38 | 252 | 815 |
| 39 | 253 | 650 |
It seems that our frequency values are traversing several orders of magnitude. Just for the sake of easy visualization, let’s leave off the top three most common values (0, 41, 2) and plot three histograms of the rest, splitting them at values below 10000, between 10000 and 12000, and above 12000.
Let’s consider these three common values for a moment, though. These values correspond to the background (i.e. unlabeled entries), left cerebral white matter, and right cerebral white matter. The extremely large value associated with “background” (15382362) is indicative of the small proportion of the data volume that is actually associated with labeled entries. Likewise, the number of voxels for left and right cerebral white matter (268665 and 265864) are indicative of the large amount of brain “real estate” occupied by the undifferentiated white matter labels. Ultimately, it will be our task as white matter segmenters to divide those areas up more sensibly, but for now let’s look at a plot of the relative size of the brain volumes occupied by these labels. We have opted not to plot the data for background, left and right cerebral white matter because of how they would dwarf the other values. Also, keep in mind, that given the resolution of the atlas we are using (1 mm^3, isotropic) the values you note in the subsequent plot can be interpreted as cubic millimeters or as milliliters (so as to compare to intuitive reference volumes like a standard 12oz soda can, which contains ~ 355 ml of liquid).
#figure out a way to plot/illustrate this better?
#drop the 0,41,2 values
labelsSubset=AtlasLabelsDataframe.loc[~AtlasLabelsDataframe.loc[:,0].isin([0,41,2])]
#reset the index
labelsSubset=labelsSubset.reset_index(drop=True)
#rename the column
labelsSubset=labelsSubset.rename(columns={0:"Label Number"})
#now split the values
lowerLabels=labelsSubset.loc[labelsSubset.iloc[:,0].le(10000)]
#reset the index
lowerLabels=lowerLabels.reset_index(drop=True)
midLabels=labelsSubset.loc[np.logical_and(~labelsSubset.iloc[:,0].le(10000),~labelsSubset.iloc[:,0].ge(12000))]
#reset the index
midLabels=midLabels.reset_index(drop=True)
upperLabels=labelsSubset.loc[labelsSubset.iloc[:,0].ge(12000)]
#reset the index
upperLabels=upperLabels.reset_index(drop=True)
#import seaborn and matplotlib for plotting
import seaborn as sns
import matplotlib.pyplot as plt
%matplotlib inline
fig, (ax1,ax2,ax3)=plt.subplots(1, 3)
fig.tight_layout(pad=3.0)
sns.countplot(y='Label Number',data=lowerLabels,ax=ax1)
#plt.gcf().set_size_inches(10,30)
#plt.figure(figsize=(10,30),dpi=200)
sns.countplot(y='Label Number',data=midLabels,ax=ax2)
#plt.gcf().set_size_inches(10,30)
#plt.figure(figsize=(10,30),dpi=200)
sns.countplot(y='Label Number',data=upperLabels,ax=ax3)
plt.gcf().set_size_inches(10,30)
#plt.figure(figsize=(10,30),dpi=200)
In the above histograms we get a sense of how many voxels (on the count axis) of the atlas correspond to each numerical label (on the Label Number’ axis). We see that some of these occur many thousands of times, while others occur quite a bit more infrequently. Furthermore, you might also notice that the bottom two plots are roughly identical, this is because these two sets of numbers (11xyz and 12xyz) correspond to the left and right variants of the same structures in the current brain. This is a standard feature of this particular freesurfer parcellation. Let’s consider the meaning of these numeric labels a bit more in depth.
Previously, when we were looking at the “paint-by-numbers” map of the United States, we noted that a parcellation is only really useful if we know the categories/framework the parcellations are mapping. As it turns out the parcellation we are looking at here is very well documented, and the meaning of these labels is quite clear. Freesurfer’s documentation includes a listing of all possible label values. Let’s load this list up here and consider the names of the regions that we have labeled here.
#this code ensures that we can navigate the WiMSE repo across multiple systems
import subprocess
import os
#get top directory path of the current git repository, under the presumption that
#the notebook was launched from within the repo directory
gitRepoPath=subprocess.check_output(['git', 'rev-parse', '--show-toplevel']).decode('ascii').strip()
#move to the top of the directory
os.chdir(gitRepoPath)
#establish path to table
FSTablePath=os.path.join(gitRepoPath,'exampleData','FreesurferLookup.csv')
#read table using pandas
FSTable=pd.read_csv(FSTablePath)
#print the table interactively
#itables.show(FSTable,columnDefs=[{"width": "10px", "targets": [0,2, 3, 4,5]}])
FSTable.tail(30)
| #No. | LabelName: | R | G | B | A | |
|---|---|---|---|---|---|---|
| 1262 | 14146 | wm_rh_S_central | 221 | 20 | 10 | 0 |
| 1263 | 14147 | wm_rh_S_cingul-Marginalis | 221 | 20 | 100 | 0 |
| 1264 | 14148 | wm_rh_S_circular_insula_ant | 221 | 60 | 140 | 0 |
| 1265 | 14149 | wm_rh_S_circular_insula_inf | 221 | 20 | 220 | 0 |
| 1266 | 14150 | wm_rh_S_circular_insula_sup | 61 | 220 | 220 | 0 |
| 1267 | 14151 | wm_rh_S_collat_transv_ant | 100 | 200 | 200 | 0 |
| 1268 | 14152 | wm_rh_S_collat_transv_post | 10 | 200 | 200 | 0 |
| 1269 | 14153 | wm_rh_S_front_inf | 221 | 220 | 20 | 0 |
| 1270 | 14154 | wm_rh_S_front_middle | 141 | 20 | 100 | 0 |
| 1271 | 14155 | wm_rh_S_front_sup | 61 | 220 | 100 | 0 |
| 1272 | 14156 | wm_rh_S_interm_prim-Jensen | 141 | 60 | 20 | 0 |
| 1273 | 14157 | wm_rh_S_intrapariet_and_P_trans | 143 | 20 | 220 | 0 |
| 1274 | 14158 | wm_rh_S_oc_middle_and_Lunatus | 101 | 60 | 220 | 0 |
| 1275 | 14159 | wm_rh_S_oc_sup_and_transversal | 21 | 20 | 140 | 0 |
| 1276 | 14160 | wm_rh_S_occipital_ant | 61 | 20 | 180 | 0 |
| 1277 | 14161 | wm_rh_S_oc-temp_lat | 221 | 140 | 20 | 0 |
| 1278 | 14162 | wm_rh_S_oc-temp_med_and_Lingual | 141 | 100 | 220 | 0 |
| 1279 | 14163 | wm_rh_S_orbital_lateral | 221 | 100 | 20 | 0 |
| 1280 | 14164 | wm_rh_S_orbital_med-olfact | 181 | 200 | 20 | 0 |
| 1281 | 14165 | wm_rh_S_orbital-H_Shaped | 101 | 20 | 20 | 0 |
| 1282 | 14166 | wm_rh_S_parieto_occipital | 101 | 100 | 180 | 0 |
| 1283 | 14167 | wm_rh_S_pericallosal | 181 | 220 | 20 | 0 |
| 1284 | 14168 | wm_rh_S_postcentral | 21 | 140 | 200 | 0 |
| 1285 | 14169 | wm_rh_S_precentral-inf-part | 21 | 20 | 240 | 0 |
| 1286 | 14170 | wm_rh_S_precentral-sup-part | 21 | 20 | 200 | 0 |
| 1287 | 14171 | wm_rh_S_suborbital | 21 | 20 | 60 | 0 |
| 1288 | 14172 | wm_rh_S_subparietal | 101 | 60 | 60 | 0 |
| 1289 | 14173 | wm_rh_S_temporal_inf | 21 | 180 | 180 | 0 |
| 1290 | 14174 | wm_rh_S_temporal_sup | 223 | 220 | 60 | 0 |
| 1291 | 14175 | wm_rh_S_temporal_transverse | 221 | 60 | 60 | 0 |
Here we see the list of possible label numbers (under #No), the the names of the anatomical regions they correspond to (under LabelName), and the internal color table assigned by Freesurfer (for plotting applications in FreeSurfer) under R,G, & B. Although there are over 1,000 entries here, any given parcellation will only have a subset of these. This is because this table (taken from the FreeSurfer website) actually covers several, mutually exclusive, parcellation schemas. In this way, two labels/numbers may correspond to overlapping areas of cortex (and thus the same constituent voxels of the associated T1 image), but not in the same data structure. For example, this would be the case with the label for the superior frontal gyrus from the Desikan-Killiany atlas[citation]–1028– the superior frontal gyrus from the “Destrieux” 2009 atlas–11116. To reiterate though, in any given parcellation/atlas data object, any particular voxel could only be associated with one of these numbers/labels, and this would be determined by the parcellation schema (i.e. “ontology”) that the parcellation data object was derived from.
Let’s limit our view of the table to consider only those labels that are present in our current NIfTI object and augment the table column with how many pixels correspond to that entry.
#create a boolean vector for the indexes which are in uniqueAtlasEntries
currentIndexesBool=FSTable['#No.'].isin(uniqueAtlasEntries)
#create new data frame with the relevant entries
currentParcellationEntries=FSTable.loc[currentIndexesBool]
#reset the indexes, whicj were disrupted by the previous operation
currentParcellationEntries=currentParcellationEntries.reset_index(drop=True)
#join the FS table subset with the voxel count
#note that we can drop the "Label Number" column of AtlasLabelCounts
#because the indexes of both dataframes track the label number
currentParcellationEntries=currentParcellationEntries.join(AtlasLabelCounts["voxel count"],how="left")
#display the table interactively
#itables.show(currentParcellationEntries,columnDefs=[{"width": "10px", "targets": [2, 3, 4,5]}])
currentParcellationEntries.tail(30)
| #No. | LabelName: | R | G | B | A | voxel count | |
|---|---|---|---|---|---|---|---|
| 162 | 12146 | ctx_rh_S_central | 221 | 20 | 10 | 0 | 4098 |
| 163 | 12147 | ctx_rh_S_cingul-Marginalis | 221 | 20 | 100 | 0 | 2303 |
| 164 | 12148 | ctx_rh_S_circular_insula_ant | 221 | 60 | 140 | 0 | 1300 |
| 165 | 12149 | ctx_rh_S_circular_insula_inf | 221 | 20 | 220 | 0 | 1997 |
| 166 | 12150 | ctx_rh_S_circular_insula_sup | 61 | 220 | 220 | 0 | 2691 |
| 167 | 12151 | ctx_rh_S_collat_transv_ant | 100 | 200 | 200 | 0 | 1942 |
| 168 | 12152 | ctx_rh_S_collat_transv_post | 10 | 200 | 200 | 0 | 665 |
| 169 | 12153 | ctx_rh_S_front_inf | 221 | 220 | 20 | 0 | 3777 |
| 170 | 12154 | ctx_rh_S_front_middle | 141 | 20 | 100 | 0 | 3542 |
| 171 | 12155 | ctx_rh_S_front_sup | 61 | 220 | 100 | 0 | 4466 |
| 172 | 12156 | ctx_rh_S_interm_prim-Jensen | 141 | 60 | 20 | 0 | 773 |
| 173 | 12157 | ctx_rh_S_intrapariet_and_P_trans | 143 | 20 | 220 | 0 | 4233 |
| 174 | 12158 | ctx_rh_S_oc_middle_and_Lunatus | 101 | 60 | 220 | 0 | 1214 |
| 175 | 12159 | ctx_rh_S_oc_sup_and_transversal | 21 | 20 | 140 | 0 | 1709 |
| 176 | 12160 | ctx_rh_S_occipital_ant | 61 | 20 | 180 | 0 | 1285 |
| 177 | 12161 | ctx_rh_S_oc-temp_lat | 221 | 140 | 20 | 0 | 2300 |
| 178 | 12162 | ctx_rh_S_oc-temp_med_and_Lingual | 141 | 100 | 220 | 0 | 3153 |
| 179 | 12163 | ctx_rh_S_orbital_lateral | 221 | 100 | 20 | 0 | 1078 |
| 180 | 12164 | ctx_rh_S_orbital_med-olfact | 181 | 200 | 20 | 0 | 1198 |
| 181 | 12165 | ctx_rh_S_orbital-H_Shaped | 101 | 20 | 20 | 0 | 2839 |
| 182 | 12166 | ctx_rh_S_parieto_occipital | 101 | 100 | 180 | 0 | 3142 |
| 183 | 12167 | ctx_rh_S_pericallosal | 181 | 220 | 20 | 0 | 1214 |
| 184 | 12168 | ctx_rh_S_postcentral | 21 | 140 | 200 | 0 | 4390 |
| 185 | 12169 | ctx_rh_S_precentral-inf-part | 21 | 20 | 240 | 0 | 2371 |
| 186 | 12170 | ctx_rh_S_precentral-sup-part | 21 | 20 | 200 | 0 | 3256 |
| 187 | 12171 | ctx_rh_S_suborbital | 21 | 20 | 60 | 0 | 629 |
| 188 | 12172 | ctx_rh_S_subparietal | 101 | 60 | 60 | 0 | 2190 |
| 189 | 12173 | ctx_rh_S_temporal_inf | 21 | 180 | 180 | 0 | 2360 |
| 190 | 12174 | ctx_rh_S_temporal_sup | 223 | 220 | 60 | 0 | 11805 |
| 191 | 12175 | ctx_rh_S_temporal_transverse | 221 | 60 | 60 | 0 | 474 |
There are several features to note about this specific parcellation schema, the “Destrieux” 2009 atlas:
The “Unknown” entries are extremely numerous. These constitute the background voxels–those for which no explicit anatomical mapping exists. Thus, even voxels without a direct mapping still have a numerical label attached to them.
Numbers below 1000 (at least of the entries we see) correspond to subcortical structures
Similarly, numbers above 10000 correspond to cortical structures
5 digit numbers that start with a 11 correspond to left hemisphere, while 5 digit numbers that start with a 12 correspond to the right hemisphere
These characteristics are not true of all parcellation schemas, but they do make a parcellation easier to work with. Particularly, numerical regularities that indicate anatomical features (i.e. cortical vs. subcortical, left hemisphere vs right hemisphere, etc.) are of great use.
Now that we have taken a look at our atlas’s labels and their voxel counts, lets see how they are spatially laid out in the data volume.
from niwidgets import NiftiWidget
#in order to have visually distinguishable areas we have to renumber the labels in the data object
#this is because niwidgets scales the color map via the min and max values of the labeling scheme,
#rather than by unique values
relabeledAtlas=atlasData.copy()
#iterate across unique label entries
for iLabels in range(len(uniqueAtlasEntries)):
#replace the current uniqueAtlasEntries value with the iLabels value
#constitutes a simple renumbering schema
relabeledAtlas[relabeledAtlas==uniqueAtlasEntries[iLabels]]=iLabels
#this code ensures that we can navigate the WiMSE repo across multiple systems
import subprocess
import os
#get top directory path of the current git repository, under the presumption that
#the notebook was launched from within the repo directory
gitRepoPath=subprocess.check_output(['git', 'rev-parse', '--show-toplevel']).decode('ascii').strip()
#move to the top of the directory
os.chdir(gitRepoPath)
#establish path to new nifti file
newAtlasPath=os.path.join(gitRepoPath,'exampleData','renumberedAtlas.nii.gz')
#store the modified atlas data in a nifti object
renumberedAtlasNifti=nib.Nifti1Image(relabeledAtlas, atlasImg.affine, atlasImg.header)
#save the object down
nib.save(renumberedAtlasNifti, newAtlasPath)
#plot it
atlas_widget = NiftiWidget(newAtlasPath)
atlas_widget.nifti_plotter(colormap='nipy_spectral')
<Figure size 432x288 with 0 Axes>
As you use the above sliders to move through the parcellation there are several things to take note of. First and perhaps most striking, while the sagittal view may be appropriately labeled, it is oriented in a somewhat disorienting fashion (typically we might expect this to be oriented such that the left and right axes of the plot correspond to the anterior-posterior axis of the brain). The coronal and axial slices are switched, but even switching them back wouldn’t be a complete solution, as the mislabeled axial slice appears to be upside down. These observations combined suggest that the affine matrix is either incorrect or is being interpreted incorrectly by the visualization function. At the moment, there’s not much we can do about this, so let’s move on to consider some other features of the visualization.
Another feature of the visualization we can note is the qualitative fashion which the parcellation actually “parcellates” (divides into pieces) the brain. The vibrant, “color-by-numbers” image we’re seeing above illustrates the distinct numerical entries in the NIfTI data object. Just as in the digital image parcellation case, each data structure element can hold a single integer value, and can thus be associated with one, and only one, label, and by extension, color. Multiple integers (and thus assignments) simply won’t fit in the same data entry in the parcellation data structure.
As a general note that is particular to brain parcellations, the anatomical entities represented in a parcellation will vary depending on what method and/or “ontology” (list of anatomical entities which are considered to “exist” for the purposes of a given parcellation) is applied. For example, the brain regions represented in the parcellation presented above are derived from the “Destrieux” 2009 atlas. The parcelled regions in this atlas (and the majority of other brain atlases, for that matter) correspond to brain regions which exhibit structural, functional, and/or microstructural homogenities relative to their neighboring brain regions (Eickhoff, Yeo, Genon, et al. 2018).
Now that we have explored three dimensional parcellations a bit, let’s move on to a consideration of these NIfTI data objects are different from tractography data objects.