Edited by: William T. Katz, Janelia Research Campus, United States
Reviewed by: Mark Kittisopikul, Janelia Research Campus, United States; Kevin Boergens, Paradromics, Inc., United States
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
Neuromorphology is crucial to identifying neuronal subtypes and understanding learning. It is also implicated in neurological disease. However, standard morphological analysis focuses on macroscopic features such as branching frequency and connectivity between regions, and often neglects the internal geometry of neurons. In this work, we treat neuron trace points as a sampling of differentiable curves and fit them with a set of branching B-splines. We designed our representation with the Frenet-Serret formulas from differential geometry in mind. The Frenet-Serret formulas completely characterize smooth curves, and involve two parameters, curvature and torsion. Our representation makes it possible to compute these parameters from neuron traces in closed form. These parameters are defined continuously along the curve, in contrast to other parameters like tortuosity which depend on start and end points. We applied our method to a dataset of cortical projection neurons traced in two mouse brains, and found that the parameters are distributed differently between primary, collateral, and terminal axon branches, thus quantifying geometric differences between different components of an axonal arbor. The results agreed in both brains, further validating our representation. The code used in this work can be readily applied to neuron traces in SWC format and is available in our open-source Python package
Not long after scientists like Ramon y Cajal started studying the nervous system with staining and microscopy, neuron morphology became a central topic in neuroscience (Parekh and Ascoli,
Currently, studying neuron morphology typically involves imaging one or more neurons, then tracing the cells and storing the traces in a digital format. Several recent initiatives have accumulated large datasets of neuron traces to facilitate morphology research. NeuroMorpho.Org, for example, hosts a total of over 140,000 neuron traces from a variety of animal species (Ascoli et al.,
Many scientists analyze neuron morphology by computing various summary features such as number of branch points, total length, and total encompassed volume. Neurolucida, a popular neuromorphology software, employs this technique. Another approach focuses on neuron topology, and uses metrics such as tree edit distance (Heumann and Wittum,
In this work, we look at neuron traces through the lens of differential geometry. In particular, we establish a system of fitting interpolating splines to the neuron traces, and computing their curvature and torsion properties. To our knowledge, curvature and torsion have never been measured in neuron traces. We applied this method to cortical projection neuron traces from two mouse brains in the MouseLight dataset from HHMI Janelia (Winnubst et al.,
First, the neuron traces were split into segments by recursively identifying the longest root to leaf path (
Next, a B-spline was fit to each point sequence using
Splines are fit to data by solving a constrained optimization problem, where a smoothing term is minimized while keeping the residual error under a specified value (Dierckx,
Sequences of fewer than 5 points, however, required lower degree splines to fully constrain the fitting procedure. For a sequence of 3 <
We recall that B-splines are not required to be parameterized by the arclength of the curve. Here, we set ξ = {0, …,
An important advantage of B-splines is that their derivatives can be computed in closed form. In fact, their derivatives are defined in terms of B-splines as shown below in Theorem 3 from Kunoth et al. (
Curvature and torsion can be easily computed because of this property. For a thrice differentiable curve
defined with the standard Euclidean norm || · ||, inner product 〈·,·〉, and cross product ×. When curvature vanishes, we define torsion to be zero as well, since the torsion equation becomes undefined. The units of curvature and torsion are both inverse length. In this work, neuron traces have units of microns, so curvature and torsion both have units of (μ
Curvature measures how much a curve deviates from being straight, and torsion measures how much a curve deviates from being planar. Together, these quantities parametrize the Frenet-Serret formulas of differential geometry. These formulas completely characterize continuously differentiable curves in three-dimensional Euclidean space, up to rigid motion (Grenander et al.,
We applied our methods to a collection of cortical projection neuron axon traces from two mouse brains in the HHMI Janelia MouseLight dataset. The precision of the reconstructions is limited by the resolution of the original two-photon block-face images, which was 0.3μ
After fitting splines to these traces, curvature and torsion magnitude were sampled every 1μ
Our first goal was to identify the length scale at which straight axon segments remain straight and curved axon segments remain curved, so we studied the autocorrelation of curvature and torsion magnitude along the axon segments. For each axon segment, the autocorrelation functions of curvature and torsion were computed along the length of the segment, yielding a collection of autocorrelation functions for each brain. Then, we evaluated whether autocorrelation at a particular lag was significantly higher than 0.3 using a one-sided
It is worth noting that, by the nature of the spline fitting procedure in Virtanen et al. (
Our second goal in the analysis was to compare curvature/torsion between segment classes. First, we estimated each segment's average curvature/torsion magnitude by taking the mean from all points that were sampled on that segment.
In order to compare different segment classes, we developed a paired sample method for testing for differences in average curvature/torsion. Different neurons represented different samples, and the average curvature/torsion of two segment classes (primary vs. collateral, collateral vs. terminal, primary vs. terminal) represented the paired measurements.
Define the random variable
We tested these hypotheses using the sign test (Neuhauser,
We also wanted to study whether these results would hold if the traces were perturbed. In particular, since the annotators vary the distance between points in their trace, we decided to randomly remove trace points and repeat the curvature/torsion measurements. Since the traces are tree structures, a trace point can be removed after connecting its child node(s) to its parent node. We produced 20 copies of the original dataset and, in each case, removed every trace point with 10% probability.
The autocorrelation functions for all segments of a brain were averaged, and they are shown in
Autocorrelation of curvature and torsion magnitude averaged across all axon segments with ±1σ confidence intervals. Curvature and torsion were sampled at every 1μ
The distributions of mean curvature and torsion are shown in
When we applied the same testing procedure to the 20 datasets with trace points randomly removed, the null hypotheses were also all rejected, in the same directions, in all cases.
The distributions of average curvature and average torsion differed between the different segment classes as shown in these kernel density estimates (which integrate to one, and therefore density has the units of μ
Neuron counts for all 36 possible curvature/torsion orderings across classes are shown in
For each neuron, average curvature and torsion was computed for all three segment classes (P, primary; C, collateral; T, terminal) and compared between classes. These heatmaps show the neuron counts for all 36 possible orderings of curvature/torsion. The most common ordering was collateral > terminal > primary for curvature and collateral > primary > terminal for torsion.
In the
Our work proposes a model of neuron morphology using continuously differentiable B-splines. From these curves, it is possible to measure kinematic properties of neuronal processes, including curvature and torsion. These techniques are freely available in our open source Python package
In most contemporary neuromorphological analysis, neuron traces are regarded as piecewise linear structures, which precludes any analysis of higher order derivatives. Our spline representation makes it possible to estimate higher order derivatives and study parameters like curvature and torsion of neuron branches. In the popular piecewise linear representation, curvature and torsion would be zero along the line segments, and undefined where the line segments meet. We simulated a piecewise linear representation by modifying our spline fitting procedure to only produce splines of degree one. Indeed, with this less sophisticated representation, curvature and torsion vanished everywhere, making them not meaningful.
Tortuosity index captures similar information to our curvature/torsion measurements and is popular in neuromorphological analysis (Stepanyants et al.,
Our methods for fitting splines and measuring curvature and torsion can be applied in neuromorphological analysis in a variety of ways, but we highlight two applications here, on a dataset of 230 projection neuron traces from two different mouse brains. We found that the autocorrelation functions of both curvature and torsion showed statistically significant correlations above 0.3 within lags of approximately 2 microns (specific lag values given in section 3.1). Next, we defined segments as either “primary,” “collateral,” or “terminal,” and found significant differences in the distributions of curvature and torsion between these classes.
The statistical analysis approach described in section 2.5 satisfies two desirable properties. First, by averaging measurements across segment classes, and pairing the data, we did not have to assume independence between segments of the same neuron. Assuming independence seemed inappropriate because, for example, segments that are connected to each other may have correlated geometry. Second, it avoided any parametric assumptions of the data, such as assuming normality of curvature/torsion measurements. A normality assumption seemed inappropriate for several reasons, including the fact that curvature is nonnegative, and that curvature/torsion was identically 0 for short segments with only 2 trace points.
Together, these findings suggest that the geometry of primary axon branches is different than that of higher order branches, such as the segments in terminal arborizations. In particular, higher order branches (collaterals and terminals) had higher curvature than primary branches. Collateral branches also had the highest torsion, but primary branches had higher torsion than terminal segments.
The primary limitation of our work is that our process of splitting a neuron trace into segments may not partition an axonal arbor into the most meaningful segment classes. This is because we needed an unambiguous classification system, while most definitions used in neuroscience literature are subjective and qualitative. For example, collaterals are broadly defined as branches that split off their parent branch at sharp angles, and arborize in a different location from other branches (Rockland,
Previous research has already indicated differences in axon geometry across neuronal cell types. For example, Stepanyants et al. (
It is also worth noting that this is not the first work to model neuron traces as continuous curves in ℝ3. For example, Duncan et al. (
This method could also be applied to measure curvature and torsion of dendrites, since dendrites also have a tree structure and are commonly stored in SWC format. However, the segment classes that we define (primary, collateral and terminal) would be inappropriate for dendrites. A segmentation classification system for dendrites would likely depend on the neuron type being studied. For example, a natural classification system of dendrites in pyramidal cells may separate apical dendrites from basal ones while dendrites in Purkinje cells would not have such a division. The researcher would have to define the dendrite segment classes according to the dataset, and the goals of the research.
It is well known that axons are pruned and modified over time (Portera-Cailliau et al.,
The datasets analyzed for this study can be found in the Open Neurodata AWS account (
MM and DT advised on the theoretical direction of the manuscript. UM advised on the data application experiments. JV advised on the presentation of the results. TA and JT designed the study, implemented the software, and managed the manuscript text/figures. All authors contributed to manuscript revision.
MM own Anatomy Works with the arrangement being managed by Johns Hopkins University in accordance with its conflict of interest policies. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
We thank the MouseLight team at HHMI Janelia for providing us with access to this data, and answering our questions about it.
The Supplementary Material for this article can be found online at:
The above plots show the relationship between segment length, and mean curvature or torsion in each segment class and brain. Each data point represents a single axon segment, and average curvature and torsion was computed by sampling the segments at a uniform spacing of 1 μ