Exploring the atlases

There are four different atlas types in ConWhat, corresponding to the 2 ontology types (Tract-based / Connectivity-Based) and 2 representation types (Volumetric / Streamlinetric).

(More on this schema here)

>>> # ConWhAt stuff
>>> from conwhat import VolConnAtlas,StreamConnAtlas,VolTractAtlas,StreamTractAtlas
>>> from conwhat.viz.volume import plot_vol_scatter

>>> # Neuroimaging stuff
>>> import nibabel as nib
>>> from nilearn.plotting import plot_stat_map,plot_surf_roi

>>> # Viz stuff
>>> %matplotlib inline
>>> from matplotlib import pyplot as plt
>>> import seaborn as sns

>>> # Generic stuff
>>> import glob, numpy as np, pandas as pd, networkx as nx

We’ll start with the scale 33 lausanne 2008 volumetric connectivity-based atlas.

Define the atlas name and top-level directory location

>>> atlas_dir = '/scratch/hpc3230/Data/conwhat_atlases'
>>> atlas_name = 'CWL2k8Sc33Vol3d100s_v01'

Initialize the atlas class

>>> vca = VolConnAtlas(atlas_dir=atlas_dir + '/' + atlas_name,
                        atlas_name=atlas_name)

loading file mapping
loading vol bbox
loading connectivity

This atlas object contains various pieces of general information

>>> vca.atlas_name

'CWL2k8Sc33Vol3d100s_v01'

>>> vca.atlas_dir

'/scratch/hpc3230/Data/conwhat_atlases/CWL2k8Sc33Vol3d100s_v01'

Information about each atlas entry is contained in the vfms attribute, which returns a pandas dataframe

Additionally, connectivity-based atlases also contain a networkx graph object vca.Gnx, which contains information about each connectome edge

>>> vca.Gnx.edges[(10,35)]

{'attr_dict': {'4dvolind': nan,
  'fullname': 'L_paracentral_to_L_caudate',
  'idx': 1637,
  'name': '10_to_35',
  'nii_file': 'vismap_grp_11-36_norm.nii.gz',
  'nii_file_id': 1637,
  'weight': 50.240000000000002,
  'xmax': 92,
  'xmin': 61,
  'ymax': 167,
  'ymin': 75,
  'zmax': 92,
  'zmin': 62}}

Individual atlas entry nifti images can be grabbed like so

>>> img = vca.get_vol_from_vfm(1637)

getting atlas entry 1637: image file /scratch/hpc3230/Data/conwhat_atlases/CWL2k8Sc33Vol3d100s_v01/vismap_grp_11-36_norm.nii.gz

>>> plot_stat_map(img)
../_images/slice_view.png

Or alternatively as a 3D scatter plot, along with the x,y,z bounding box

>>> vca.bbox.ix[1637]

xmin     61
xmax     92
ymin     75
ymax    167
zmin     62
zmax     92
Name: 1637, dtype: int64

>>> ax = plot_vol_scatter(vca.get_vol_from_vfm(1),c='r',bg_img='nilearn_destrieux',
>>>                         bg_params={'s': 0.1, 'c':'k'},figsize=(20, 15))
>>> ax.set_xlim([0,200]); ax.set_ylim([0,200]); ax.set_zlim([0,200]);

getting atlas entry 1: image file /scratch/hpc3230/Data/conwhat_atlases/CWL2k8Sc33Vol3d100s_v01/vismap_grp_39-56_norm.nii.gz
../_images/scatter_view.png

We can also view the weights matrix like so:

>>> fig, ax = plt.subplots(figsize=(16,12))
>>> sns.heatmap(np.log1p(vca.weights),xticklabels=vca.region_labels,
>>>               yticklabels=vca.region_labels,ax=ax);
>>> plt.tight_layout()
../_images/weights_matrix.png

The vca object also contains x,y,z bounding boxes for each structure

We also stored additional useful information about the ROIs in the associated parcellation, including cortical/subcortical labels

>>> vca.cortex

array([ 1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.])

…hemisphere labels

>>> vca.hemispheres

array([ 1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,
        1.,  1.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.])

…and region mappings to freesurfer’s fsaverage brain

>>> vca.region_mapping_fsav_lh

array([ 24.,  29.,  28., ...,  16.,   7.,   7.])

>>> vca.region_mapping_fsav_rh

array([ 24.,  29.,  22., ...,   9.,   9.,   9.])

which can be used for, e.g. plotting ROI data on a surface

>>> f = '/opt/freesurfer/freesurfer/subjects/fsaverage/surf/lh.inflated'
>>> vtx,tri = nib.freesurfer.read_geometry(f)
>>> plot_surf_roi([vtx,tri],vca.region_mapping_fsav_lh);
../_images/rois_on_surf.png