Introduce a new method to calculate partition asymmetry by Uylings (#930

)

CHANGELOG.rst

@@ -1,6 +1,11 @@
 Changelog
 =========
 
+Version 2.3.0
+-------------
+- Introduce a new method to calculate partition asymmetry by Uylings. See docstring of
+  :func:`neurom.features.neuritefunc.partition_asymmetries`.
+
 Version 2.2.1
 -------------
 - Fix 'section_path_lengths' feature for Population

neurom/features/bifurcationfunc.py

@@ -99,10 +99,12 @@ def bifurcation_partition(bif_point):
     return max(n, m) / min(n, m)
 
 
-def partition_asymmetry(bif_point):
+def partition_asymmetry(bif_point, uylings=False):
     """Calculate the partition asymmetry at a bifurcation point.
 
-    Partition asymmetry is defined in https://www.ncbi.nlm.nih.gov/pubmed/18568015
+    By default partition asymmetry is defined as in https://www.ncbi.nlm.nih.gov/pubmed/18568015.
+    However if ``uylings=True`` is set then
+    https://jvanpelt.nl/papers/Uylings_Network_13_2002_397-414.pdf is used.
 
     The number of nodes in each child tree is counted. The partition
     is defined as the ratio of the absolute difference and the sum
@@ -113,9 +115,13 @@ def partition_asymmetry(bif_point):
 
     n = float(sum(1 for _ in bif_point.children[0].ipreorder()))
     m = float(sum(1 for _ in bif_point.children[1].ipreorder()))
-    if n == m:
-        return 0.0
-    return abs(n - m) / abs(n + m)
+    c = 0
+    if uylings:
+        c = 2
+        if n + m <= c:
+            raise NeuroMError('Partition asymmetry cant be calculated by Uylings because the sum of'
+                              'terminal tips is less than 2.')
+    return abs(n - m) / abs(n + m - c)
 
 
 def partition_pair(bif_point):

neurom/features/neuritefunc.py

@@ -372,28 +372,27 @@ def remote_bifurcation_angles(neurites, neurite_type=NeuriteType.all):
                          iterator_type=Section.ibifurcation_point)
 
 
-@feature(shape=(...,), name='partition')
-def bifurcation_partitions(neurites, neurite_type=NeuriteType.all):
-    """Partition at bifurcation points of a collection of neurites."""
-    return map(bifurcationfunc.bifurcation_partition,
-               iter_sections(neurites,
-                             iterator_type=Section.ibifurcation_point,
-                             neurite_filter=is_type(neurite_type)))
-
-
 @feature(shape=(...,), name='partition_asymmetry')
-def partition_asymmetries(neurites, neurite_type=NeuriteType.all, variant='branch-order'):
+def partition_asymmetries(neurites,
+                          neurite_type=NeuriteType.all,
+                          variant='branch-order',
+                          method='petilla'):
     """Partition asymmetry at bifurcation points of a collection of neurites.
 
     Variant: length is a different definition, as the absolute difference in
     downstream path lenghts, relative to the total neurite path length
+    Method: 'petilla' or 'uylings'. The former is default. The latter uses ``-2`` shift. See
+    :func:`neurom.features.bifurcationfunc.partition_asymmetry`
     """
     if variant not in {'branch-order', 'length'}:
-        raise ValueError('Please provide a valid variant for partition asymmetry,\
-                         found %s' % variant)
+        raise ValueError('Please provide a valid variant for partition asymmetry,'
+                         f'found {variant}')
+    if method not in {'petilla', 'uylings'}:
+        raise ValueError('Please provide a valid method for partition asymmetry,'
+                         'either "petilla" or "uylings"')
 
     if variant == 'branch-order':
-        return map(bifurcationfunc.partition_asymmetry,
+        return map(partial(bifurcationfunc.partition_asymmetry, uylings=method == 'uylings'),
                    iter_sections(neurites,
                                  iterator_type=Section.ibifurcation_point,
                                  neurite_filter=is_type(neurite_type)))
@@ -410,6 +409,15 @@ def partition_asymmetries(neurites, neurite_type=NeuriteType.all, variant='branc
     return asymmetries
 
 
+@feature(shape=(...,), name='partition')
+def bifurcation_partitions(neurites, neurite_type=NeuriteType.all):
+    """Partition at bifurcation points of a collection of neurites."""
+    return map(bifurcationfunc.bifurcation_partition,
+               iter_sections(neurites,
+                             iterator_type=Section.ibifurcation_point,
+                             neurite_filter=is_type(neurite_type)))
+
+
 # Register `partition_asymmetries` variant
 _partition_asymmetry_length = partial(partition_asymmetries, variant='length')
 update_wrapper(_partition_asymmetry_length, partition_asymmetries)  # this fixes the docstring

tests/features/test_bifurcationfunc.py

@@ -37,9 +37,10 @@
 import neurom as nm
 from neurom import load_neuron
 from neurom.exceptions import NeuroMError
-
 from neurom.features import bifurcationfunc as bf
 
+import pytest
+
 DATA_PATH = Path(__file__).parent.parent / 'data'
 SWC_PATH = DATA_PATH / 'swc'
 SIMPLE = nm.load_neuron(SWC_PATH / 'simple.swc')
@@ -78,6 +79,9 @@ def test_partition_asymmetry():
     root = SIMPLE2.neurites[0].root_node
     assert bf.partition_asymmetry(root) == 0.5
     assert bf.partition_asymmetry(root.children[0]) == 0.0
+    assert bf.partition_asymmetry(root, True) == 1.0
+    with pytest.raises(NeuroMError, match='Uylings'):
+        bf.partition_asymmetry(root.children[0], True)
 
     leaf = root.children[0].children[0]
     assert_raises(NeuroMError, bf.partition_asymmetry, leaf)

tests/features/test_neuritefunc.py

@@ -313,7 +313,10 @@ def test_partition_asymmetry():
                     (0.0625, 0.06666666666666667))
 
     with pytest.raises(ValueError):
-        _nf.partition_asymmetries(SIMPLE, variant='unvalid-variant')
+        _nf.partition_asymmetries(SIMPLE, variant='invalid-variant')
+
+    with pytest.raises(ValueError):
+        _nf.partition_asymmetries(SIMPLE, method='invalid-method')
 
 
 def test_segment_lengths():