Refactor and fix peak_local_max (#4760)

* Remove relabelling

* Add trivial case support in label regions

* Fix min_distance support

* Fix test_peak

* Fix test_corner

* Fix corner_peaks doctest

* Vectorize KDtree query

* Update skimage/feature/peak.py

Co-authored-by: Marianne Corvellec <marianne.corvellec@ens-lyon.org>

* Fix merge conflicts

* Fix docstring

* Refactor _get_high_intensity_peaks

* Add ensure_spacing function

* Update _get_threshold

* Update _exclude_border

* Remove warning in test

* Add docstring for _get_threshold

* Add some comments to ease review

* Fix test_input_labels_unmodified

* Add note about behaviour change in 0.18

* Fix ensure_spacing

* Revert test_peak.test_empty

* Fix ensure_spacing

* Revert corner_peaks doctest

* Revert test_disk

* Simplify test_input_labels_unmodified

* Revert warning msg when footprint.size < 2

* Update skimage/_shared/tests/test_coord.py

Co-authored-by: Gregory R. Lee <grlee77@gmail.com>

* Add release notes

* Remove warning if footprint.size=1: the warning assertion stands only if min_distance < 1

Co-authored-by: Marianne Corvellec <marianne.corvellec@ens-lyon.org>
Co-authored-by: Riadh <r.fezzani@vitadx.com>
Co-authored-by: Juan Nunez-Iglesias <juan.nunez-iglesias@monash.edu>
Co-authored-by: Gregory R. Lee <grlee77@gmail.com>

doc/release/release_dev.rst

@@ -73,7 +73,7 @@ Improvements
 - ``measure.label`` has been accelerated for boolean input images, by using
   ``scipy.ndimage``'s implementation for this case (#4945).
 - ``util.apply_parallel`` now works with multichannel data (#4927).
-
+- ``skimage.feature.peak_local_max`` supports now any Minkowski distance.
 
 
 API Changes
@@ -93,6 +93,8 @@ API Changes
   descriptive name.
 - The ``level`` parameter of ``measure.find_contours`` is now a keyword
   argument, with a default value set to (max(image) - min(image)) / 2.
+- ``p_norm`` argument was added to ``skimage.feature.peak_local_max``
+  to add support for Minkowski distances.
 
 
 Bugfixes
@@ -122,6 +124,12 @@ Bugfixes
   NaNs are found in the computation (as a result of NaNs in input images).
   Before this fix, an incorrect value could be returned where the input images
   had NaNs (#4886).
+- ``min_distance`` is now enforced for ``skimage.feature.peak_local_max``
+  (#2592).
+- Peak detection in labels is fixed in ``skimage.feature.peak_local_max``
+  (#4756).
+- Input ``labels`` argument renumbering in ``skimage.feature.peak_local_max``
+  is avoided (#5047).
 
 
 Deprecations

skimage/_shared/coord.py

@@ -0,0 +1,51 @@
+import numpy as np
+from scipy.spatial import cKDTree, distance
+
+
+def ensure_spacing(coord, spacing=1, p_norm=np.inf):
+    """Returns a subset of coord where a minimum spacing is guaranteed.
+
+    Parameters
+    ----------
+    coord : ndarray
+        The coordinates of the considered points.
+    spacing : float
+        the maximum allowed spacing between the points.
+    p_norm : float
+        Which Minkowski p-norm to use. Should be in the range [1, inf].
+        A finite large p may cause a ValueError if overflow can occur.
+        ``inf`` corresponds to the Chebyshev distance and 2 to the
+        Euclidean distance.
+
+    Returns
+    -------
+    output : ndarray
+        A subset of coord where a minimum spacing is guaranteed.
+
+    """
+
+    output = coord
+    if len(coord):
+        # Use KDtree to find the peaks that are too close to each other
+        tree = cKDTree(coord)
+
+        indices = tree.query_ball_point(coord, r=spacing, p=p_norm)
+        rejected_peaks_indices = set()
+        for idx, candidates in enumerate(indices):
+            if idx not in rejected_peaks_indices:
+                # keep current point and the points at exactly spacing from it
+                candidates.remove(idx)
+                dist = distance.cdist([coord[idx]],
+                                      coord[candidates],
+                                      distance.minkowski,
+                                      p=p_norm).reshape(-1)
+                candidates = [c for c, d in zip(candidates, dist)
+                              if d < spacing]
+
+                # candidates.remove(keep)
+                rejected_peaks_indices.update(candidates)
+
+        # Remove the peaks that are too close to each other
+        output = np.delete(coord, tuple(rejected_peaks_indices), axis=0)
+
+    return output

skimage/_shared/tests/test_coord.py

@@ -0,0 +1,43 @@
+import pytest
+import numpy as np
+from scipy.spatial.distance import pdist, minkowski
+from skimage._shared.coord import ensure_spacing
+
+
+@pytest.mark.parametrize("p", [1, 2, np.inf])
+def test_ensure_spacing_trivial(p):
+    # --- Empty input
+    assert ensure_spacing([], p_norm=p) == []
+
+    # --- Verified spacing
+    coord = np.random.randn(100, 2)
+
+    # --- 0 spacing
+    assert np.array_equal(coord, ensure_spacing(coord, spacing=0, p_norm=p))
+
+    # Spacing is chosen to be half the minimum distance
+    spacing = pdist(coord, minkowski, p=p).min() * 0.5
+
+    out = ensure_spacing(coord, spacing=spacing, p_norm=p)
+
+    assert np.array_equal(coord, out)
+
+
+@pytest.mark.parametrize("ndim", [1, 2, 3, 4, 5])
+def test_ensure_spacing_nD(ndim):
+    coord = np.ones((5, ndim))
+
+    expected = np.ones((1, ndim))
+
+    assert np.array_equal(ensure_spacing(coord), expected)
+
+
+@pytest.mark.parametrize("p", [1, 2, np.inf])
+def test_ensure_spacing_p_norm(p):
+    coord = np.random.randn(100, 2)
+
+    # --- Consider the average distance btween the point as spacing
+    spacing = np.median(pdist(coord, minkowski, p))
+    out = ensure_spacing(coord, spacing=spacing, p_norm=p)
+
+    assert pdist(out, minkowski, p).min() > spacing