Feature map / return_indices

[kevinjohncutler]

I am curious if it would be possible to include the return_indices argument/ output of scipy.ndimage.distance_transform_edt, which returns a map of the closest object index for each pixel of the distance field. This is only relevant for the distance transform of the background, where it can be important to know which labelled object is closest.

[william-silversmith]

Hi Kevin,

I'll have to think about it. The existing algorithm doesn't propagate this information so some rethinking will be necessary. This is pretty similar to (maybe the same as) computing the voronoi diagram of a given object. Another github user requested this functionality a little while ago for doing fast morphological dilations.

Very happy you're getting a lot of mileage out of this library, I've recommended cellpose a at least once.

[acycliq]

Hi William!

Thanks for this. I just want to say that I also find the return_indices output value quite useful.

[william-silversmith]

Off the top of my head, I can think of an okay, but not great, algorithm for computing them from the output of EDT that would be relatively easy to implement. Something less than quadratic at least. I'm pretty busy at the moment, so it might take a substantial amount of time before I get a crack at it.

If you can submit a PR or a design, that might make things go faster.

Thanks so much for using EDT and being interested in its development!

[william-silversmith]

So I serendipitously discovered that the feature transform is a separable part of scipy. I might be able to include it as an optional feature if you also have scipy installed. The implementation would look something like this:

def feature_transform(labels, anisotropy):
    import scipy.ndimage._nd_image
    labels = np.atleast_1d(np.where(labels, 1, 0).astype(np.int8))
    ft = np.zeros([labels.ndim] + list(labels.shape), dtype=np.int32)
    scipy.ndimage._nd_image.euclidean_feature_transform(labels, anisotropy, ft)
    return ft

[william-silversmith]

Actually, even better you can calculate the feature transform in scipy without the distance transform.

ft = scipy.ndimage.distance_transform_edt(input, sampling, return_distances=False, return_indices=True) https://docs.scipy.org/doc/scipy/reference/generated/scipy.ndimage.distance_transform_edt.html