Returning specific numpy label indices for an image

Question

I'd like to segment my image using numpy's label and then based on the number of indices found in each label remove those which satisfy my criteria. For example if an image with regions in it that I'd segmented were created like this and segmented using scipy's label:

from numpy import ones, zeros
from numpy.random import random_integers
from scipy.ndimage import label

image = zeros((512, 512), dtype='int')
regionator = ones((11, 11), dtype='int')
xs = random_integers(5, 506, size=500)
ys = random_integers(5, 506, size=500)

for x, y in zip(xs, ys):
    image[x-5:x+6, y-5:y+6] = regionator

labels, n_labels = label(image)

Now I'd like to retrieve the indices for each region which has a size greater than 121 pixels (or one regionator size). I'd then like to take those indices and set them to zero so they are no longer part of the labeled image. What is the most efficient way to accomplish this task?

Essentially something similar to MATLAB's regionprops or utilizing IDL's reverse_indices output from its histogram function.

Answer 1

I would use bincount and threshold the result to make a lookup table:

import numpy as np

threshold = 121

size = np.bincount(labels.ravel())
keep_labels = size <= threshold
# Make sure the background is left as 0/False
keep_labels[0] = 0
filtered_labels = keep_labels[labels]

On the last above I index the array keep_labels with the array labels. This is called advanced indexing in numpy and it requires that labels be an integer array. Numpy then uses the elements of labels as indices to keep_labels and produces an array the same shape as labels.

Answer 2

Here's what I've found to work for me so far with good performance even for large datasets.

Using the get indices process taken from here I've come to this:

from numpy import argsort, histogram, reshape, where
import bisect

h = histogram(labels, bins=n_labels)
h_inds = where(h[0] > 121)[0]

labels_f = labels.flatten()
sortedind = argsort(labels_f)
sorted_labels_f = labels_f[sortedind]
inds = []
for i in range(1, len(h_inds)):
    i1 = bisect.bisect_left(sorted_labels_f, h[1][h_inds[i]])
    i2 = bisect.bisect_right(sorted_labels_f, h[1][h_inds[i]])
    inds.extend(sortedind[i1:i2])

# Now get rid of all of those indices that were part of a label
# larger than 121 pixels
labels_f[inds] = 0
filtered_labels = reshape(labels_f, (512, 512))