T02: Hemisphere Analysis#

This notebook exemplifies how to use SNT to quantify ipsilateral/contralateral projection patterns of complete neurons in the adult mouse brain (MouseLight data).

Hint

By the end of this notebook, you will be able to:

  • Load and visualize neuronal reconstructions from the MouseLight database

  • Load Allen Mouse Brain CCF data, including 3D meshes delineating brain areas

  • Perform hemisphere-specific analysis to quantify ipsilateral vs. contralateral projections

  • Use AllenUtils to work with brain atlas annotations and anatomical coordinates

  • Generate and analyze histograms of morphological measurements with statistical annotations

  • Visualize 3D neuronal morphology with anatomical context using Viewer3D

Estimated Time: 30 minutes

Note

Make sure to read these resources before running this notebook:

Setup and Initialization#

As always, we start by initializing pysnt. We’ll also set some generic options suitable for notebook execution:

import pysnt
pysnt.set_option('java.logging.level', 'Error') # Set logging level to 'Error' to reduce console output (see Overview for details)
pysnt.set_option('display.chart_format', 'svg') # SVG plots and histograms
pysnt.set_option('display.chart_dpi', 150) # Rendering resolution

import pysnt
pysnt.initialize()

Now we can import the modules that we will be needing in this notebook. Because we will be working with the Allen Mouse Brain CCF, we’ll be using AllenUtils. Note that other species/atlas have an equivalent counterpart:

from pysnt import Tree # as per tutorial 01
from pysnt.analysis import TreeStatistics # as per tutorial 01
from pysnt.annotation import AllenUtils
from pysnt.io import MouseLightLoader

Load data#

Let’s retrieve a cell from the primary motor cortex (MOp) and extract its axon:

loader = MouseLightLoader("AA0876")
if not loader.isDatabaseAvailable():
    print("Could not connect to ML database", "Error")
if not loader.idExists():
    print("Somehow the specified id was not found", "Error")
axon = loader.getTree('axon')

Let’s display it. We will assemble a 3D scene, loading some of its surrounding structures, as well as the brain contour:

from pysnt.viewer import Viewer3D

viewer = Viewer3D() # Initialize the viewer

# Add reconstruction: distinguish axon from dendrites
axon.setColor("blue")
dend = loader.getTree('dendrite')
dend.setColor("red")
viewer.add([axon, dend])

# Add meshes associated with cell
root_compartment = axon.getRoot().getAnnotation() # get brain annotation at the soma (axon root)
root_mesh = root_compartment.getMesh() # get mesh of annotation
viewer.add(root_mesh)

# Let's add some other brain regions for context
th = AllenUtils.getCompartment("Thalamus")
viewer.add(th.getMesh())

# Load brain contour, toggle light mode, and display static scene
brain_mesh = viewer.loadRefBrain('ccf')
viewer.setEnableDarkMode(False)
orthoviews = pysnt.display(viewer, orthoview=True, figsize=(18, 6))
../_images/bd8b5fa2d8e75507dc8c6b626c5a3d33da356cd0ec0fc4796a37c20b2b99a6fb.png

Let’s label the views with anatomical planes. Once more, we can use AllenUtils for this. We’ll create a function, for future re-use:

def display_orthoviews(viewer):
    from pysnt.annotation import AllenUtils


    # Retrieve orthoview layout from pysnt.display()
    orthoviews = pysnt.display(viewer, orthoview=True, figsize=(18, 6)) # SNTObject dictionary
    fig = orthoviews['data'] # access the matplotlib figure

    # Replace the titles
    for axis in fig.get_axes():
        anatomical_plane = AllenUtils.getAnatomicalPlane(axis.get_title()) # Maps "xy" to "sagittal", "yz" to "coronal", ...
        if anatomical_plane:
            axis.set_title(anatomical_plane, fontsize=14, fontweight='regular')
        # Remove the box around each view
        for spine in axis.spines.values():
            spine.set_visible(False)

    # Force redraw
    fig.canvas.draw()
    return fig


fig = display_orthoviews(viewer)
../_images/2831cbe9dc6cfed0fae39813f2c4ca000b3ee0cc634de732938f2dcdf7a36d81.png

Does the axon have contralateral projections?#

Let’s answer this programatically. First we need to retrieve the plane defining the mid-sagittal midline:

midline_axis = AllenUtils.getAxisDefiningSagittalPlane() # the cartesian axis: 0 -> X; 1 -> Y; 2 -> Z
print(f"The midline axis for the Allen CCF is: {midline_axis} ({'XYZ'[midline_axis]} axis)")
The midline axis for the Allen CCF is: 2 (Z axis)

Thus, in the Allen CCF, medio-lateral positions are mapped to the Z axis. We can confirm this by visualizing the midline plane in the scene we created before:

# retrieve the bounding box of the brain mesh
bbox = brain_mesh.getBoundingBox('both') # either 'left', 'right', or 'both' hemisphere(s)
annotated_plane = viewer.annotateMidPlane(bbox, midline_axis, "Midplane")
annotated_plane.setColor("yellow") 
fig = display_orthoviews(viewer)
../_images/a6b3f60bd5d91aaffe47db6567e7d363db0a66b637676f94de25ce73f87c7a7c.png

So, in practice we can get immediate insights into ipsilateral/contralateral hemisphere patterns by looking at the distribution of Z coordinates of reconstruction nodes:

tree_stats = TreeStatistics(axon)
hist = tree_stats.getNodeStatistics().getHistogram("z coordinates")
hist.annotateXline(AllenUtils.brainCenter().getZ(), "Midline", 'red')  # annotate Z coordinate of midline in the histogram
hist.annotateXline(axon.getRoot().getZ(), "Soma", 'blue') # annotate Z coordinate (medio-lateral position) of soma in the histogram
pysnt.display(hist)
[Run$_main] INFO smile.stat.distribution.GaussianMixture - The BIC of Mixture(1)[1.00 x Gaussian Distribution(6290.2334, 2092.8720)] = -48162.5853
[Run$_main] INFO smile.stat.distribution.GaussianMixture - The BIC of Mixture(2)[0.49 x Gaussian Distribution(5525.7847, 1964.2185) + 0.51 x Gaussian Distribution(7014.3831, 1948.0715)] = -48165.2706
../_images/565f95d914dda56123ff30c8abb371653e14da16391705a993269d5f391c85a9.png
{'type': matplotlib.figure.Figure,
 'data': <Figure size 400x400 with 1 Axes>,
 'metadata': {'source_type': 'SNTChart',
  'format': 'svg',
  'scale': 1.0,
  'is_combined': False,
  'title': 'Hist. Z coordinates',
  'containsValidData': True,
  'isLegendVisible': False},
 'error': None}

The plot shows a bimodal distribution (a distribution with two distinct peaks or modes). The first peak (mode) is associated with the left hemisphere, while the second peak (mode) is associated with the right hemisphere. The separation valley between the two peaks reflects the absence of branching near the midline.

This is corroborated by looking at the position of end-points (tips) and branch-points (junctions):

histograms = []
selected_nodes = ["Branch Points", "End Points"]
for points in selected_nodes:
    hist = tree_stats.getNodeStatistics(points).getHistogram("z coordinates")
    hist.annotateXline(AllenUtils.brainCenter().getZ(), "Midline", 'red')
    hist.annotateXline(axon.getRoot().getZ(), "Soma", 'blue')
    histograms.append(hist)

pysnt.display(histograms, panel_titles=selected_nodes)
[Run$_main] INFO smile.stat.distribution.GaussianMixture - The BIC of Mixture(1)[1.00 x Gaussian Distribution(6584.1750, 1956.9267)] = -2074.4940
[Run$_main] INFO smile.stat.distribution.GaussianMixture - The BIC of Mixture(2)[0.49 x Gaussian Distribution(5853.2888, 1884.9454) + 0.51 x Gaussian Distribution(7274.0006, 1754.4687)] = -2081.7621
[Run$_main] INFO smile.stat.distribution.GaussianMixture - The BIC of Mixture(1)[1.00 x Gaussian Distribution(6592.1846, 2048.7837)] = -2094.0926
[Run$_main] INFO smile.stat.distribution.GaussianMixture - The BIC of Mixture(2)[0.46 x Gaussian Distribution(5822.9659, 1997.5865) + 0.54 x Gaussian Distribution(7251.8523, 1843.5914)] = -2101.2720
{'type': matplotlib.figure.Figure,
 'data': <Figure size 700x350 with 2 Axes>,
 'metadata': {'source_type': 'SNTChart_List',
  'sntchart_count': 2,
  'displayed_count': 2,
  'panel_layout': 'auto',
  'title': None},
 'error': None}
../_images/bc7b4d28d56cca429e360788b3b1d1a234ceb06c9951c503e674eb6689c9446e.png

Since we are not using the 3D viewer in interactive mode, it is useful to visualize the axon by overlaying the distance of each node to the soma:

from pysnt.analysis import TreeColorMapper

temp_viewer = Viewer3D()
color_mapper = TreeColorMapper()
color_mapper.map(axon, 'distance to soma', 'plasma') # tree, metric, lut/colormap (case-sensitive)
temp_viewer.add(axon)
temp_viewer.loadRefBrain('ccf')
temp_viewer.setEnableDarkMode(False)
fig = display_orthoviews(temp_viewer)
../_images/d27aec690848aa8f754e94d055fb53e45c5f864ebba374453fe467a0f2573447.png

Does the axon target the same brain area evenly across hemispheres?#

To answer this it is convenient to first tag all the nodes in the reconstruction with ‘left’, ‘right’ hemisphere tags.

AllenUtils.assignHemisphereTags(axon) # tag each node with 'left'/'right' labels

Now when we query TreeStatistics to obtain annotated lengths (i.e., the amount of axonal cable associated with the brain areas innervated by the axon), we can split results by hemisphere. To simplify things, we are going to parse only brain areas of “mid-level ontology” (i.e., those at mid-depth in the CCF ontology, which, as of this writing, has a maximum depth of 10):

# define a function that highlights the soma compartment in a histogram
def annotate_soma(histogram, ontology_depth):
    global loader, depth
    soma_compartment = loader.getSomaCompartment()
    if soma_compartment.getOntologyDepth() > ontology_depth:
        ancestor_depth = ontology_depth - soma_compartment.getOntologyDepth()
        soma_compartment = soma_compartment.getAncestor(ancestor_depth)
    histogram.annotateCategory(soma_compartment.acronym(), "soma")

depth = 6
histograms = []
for hemisphere in ["left", "right"]:
    hist = tree_stats.getAnnotatedLengthHistogram(depth, hemisphere)
    annotate_soma(hist, depth)
    histograms.append(hist)

pysnt.display(histograms, panel_titles=["Left Hemisphere", "Right Hemisphere"])
{'type': matplotlib.figure.Figure,
 'data': <Figure size 700x350 with 2 Axes>,
 'metadata': {'source_type': 'SNTChart_List',
  'sntchart_count': 2,
  'displayed_count': 2,
  'panel_layout': 'auto',
  'title': None},
 'error': None}
../_images/5f7e89c97334114e35c593bd9c862019863c13869beadda43af093b2dcd3f6f3.png

This makes it easier to look for biases in the distributions between hemispheres. But looking at the data in this format can be rather cumbersome. It would be better to visualize both distributions side-by-side:

hist = tree_stats.getAnnotatedLengthHistogram(depth, "ratio")
hist.setFontSize(13)
annotate_soma(hist, depth)
pysnt.display(hist)
../_images/bc00ea38a40eb542d4240a7748e5777e930e9212b988afe50bcc7852324b6916.png
{'type': matplotlib.figure.Figure,
 'data': <Figure size 400x400 with 1 Axes>,
 'metadata': {'source_type': 'SNTChart',
  'format': 'svg',
  'scale': 1.0,
  'is_combined': False,
  'title': 'AA0876 (axon) Annotated Length',
  'containsValidData': True,
  'isLegendVisible': True},
 'error': None}

Now it is easier to notice that in the Isocortex

Exclusively ipsilateral projections: to the Agranular insular area (AI), Ectorhinal area (ECT)

Exclusively contralateral projections: to the Frontal Pole (FP), and Prelimbic cortex (PL)

Let’s visualize some of these regions:

brain_areas = ["ECT", "PL", "FP"]
for brain_area in brain_areas:
    mesh = AllenUtils.getCompartment(brain_area).getMesh()
    viewer.add(mesh)
viewer.removeAnnotation(annotated_plane) # remove midline plane
fig = display_orthoviews(viewer)
../_images/509830898afa1e271f9b95f7dbd8ae4588e9b32250b6c3499e60cfd8b2ad3a83.png

Followup Questions#

  1. Which kind of anatomical structures are associated with axonal processes in the ipsilateral ECT? unbranched processes? end-points? branch-points?

  2. How can we visualize the axon around that hemi-area?

Data Sources and References#

Data used in this notebook:

  • Cell AA0876 of the MouseLight database, under a Creative Commons Attribution 4.0 International License (CC BY 4.0)

  • The Allen Mouse Brain Common Coordinate Framework (CCFv3) is openly accessible at https://atlas.brain-map.org/ under the Allen Institute’s Terms of Use

See SNT citation for details on how to properly cite SNT.