Merge pull request #3 from NeuroDataDesign/compute_overlaps

remove labeling of ground truth

source/bstadt/Quality/Quality.py

@@ -1,17 +1,56 @@
 import numpy as np
 from skimage.measure import label
 
+def bounding_box(img):
+    """
+    Returns the z, y, x vectors that create a bounding box
+    of a mask.
+    """
+    z = np.any(img, axis=(1, 2))
+    y = np.any(img, axis=(0, 2))
+    x = np.any(img, axis=(0, 1))
+
+    zmin, zmax = np.where(z)[0][[0, -1]]
+    ymin, ymax = np.where(y)[0][[0, -1]]
+    xmin, xmax = np.where(x)[0][[0, -1]]
+
+    return (zmin, zmax + 1), (ymin, ymax + 1), (xmin, xmax + 1)
+
+
+def get_uniques(ar):
+    """
+    Returns an ordered numpy array of unique integers in an array.
+    This runs about four times faster than numpy.unique().
+
+    Parameters
+    ----------
+    ar : array_like
+        Input array. This will be flattened.
+
+    Returns
+    -------
+    uniques : ndarray
+        The sorted unique values.
+    """
+    bins = np.zeros(np.max(ar) + 1, dtype=int)
+    bins[ar.ravel()] = 1
+    uniques = np.nonzero(bins)[0]
+
+    return uniques
+
+
 def get_unique_overlap(foreground, background, i):
     '''
     Calculates the number of unique background labels in the foreground at i
     Does not count background label of 0
     '''
+    z, y, x = bounding_box(foreground == i)
+
+    foreground = foreground[z[0]:z[1], y[0]:y[1], x[0]:x[1]]
+    background = background[z[0]:z[1], y[0]:y[1], x[0]:x[1]]
 
-    #This runs about 4 times faster than np.unique()
     overlaps = np.multiply((foreground == i), background)
-    bins = np.zeros(np.max(overlaps) + 1, dtype=int)
-    bins[overlaps.ravel()] = 1
-    uniques = np.nonzero(bins)[0]
+    uniques = get_uniques(overlaps)
 
     num_unique = len(uniques)
 
@@ -24,24 +63,27 @@ def get_unique_overlap(foreground, background, i):
     return num_unique
 
 
-def compute_overlap_array(predictions, gt):
-
-    predictionLabels = label(predictions)
-    maxPredictionLabel = np.max(predictionLabels)
+def compute_overlap_array(predictions, gt, compare_annotations=False):
+    """
+    When comparing two annotation volumes, set compare_annotations to True.
+    """
+    if not compare_annotations:
+        predictions = label(predictions)
+        prediction_uniques = get_uniques(predictions)[1:]
+    elif compare_annotations:
+        prediction_uniques = get_uniques(predictions)[1:]
 
-    gtLabels = label(gt)
-    maxGtLabel = np.max(gtLabels)
+    gt_uniques = get_uniques(gt)[1:]
 
     #first, look at how many unique predictions
     #overlap with a single gt synapse
-    predictionPerGt = [get_unique_overlap(gtLabels, predictionLabels, i)\
-                       for i in range(1, maxGtLabel + 1)]
-
+    predictionPerGt = [get_unique_overlap(gt, predictions, i)
+                       for i in gt_uniques]
 
     #next, look at how many unique synapses overlap
     #with a single synapse prediction
-    gtPerPrediction = [get_unique_overlap(predictionLabels, gtLabels, i)\
-                       for i in range(1, maxPredictionLabel + 1)]
+    gtPerPrediction = [get_unique_overlap(predictions, gt, i)
+                       for i in prediction_uniques]
 
     return {'predictionPerGt': predictionPerGt,
             'gtPerPrediction': gtPerPrediction}

source/bstadt/Util/Util.py

@@ -1,3 +1,6 @@
+import numpy as np
+from random import randrange as rand
+
 def delta_epsilon(a, b, e):
     return abs(a-b) < e