Collectives™ on Stack Overflow

Find centralized, trusted content and collaborate around the technologies you use most.

Teams

Q&A for work

Connect and share knowledge within a single location that is structured and easy to search.

After reading this post , and play also with SciKit-image I found a difference in Python compared to MATLAB's function imregionalmax .

I have these lines of code:

from skimage.feature import peak_local_max
manos = np.ones([5,5])
manos[2,2] = 0.
manos[2,4] = 2.
giannis = peak_local_max(manos,min_distance=1, indices=False, exclude_border=False)
giorgos = ndimage.filters.maximum_filter(manos, footprint=np.ones([3,3]))
giorgos = (giorgos == manos)
I would expect a 2D array with only one True value ([2,4]) for variables giannis or giorgos as I get in MATLAB. Instead I take more than one maximum. 
Any idea why this works this way and how to make it work like in MATLAB?
Both giannis and giorgos are similar in that they find pixels that are equal or larger than the other pixels in the 3x3 neighborhood. I believe giannis would have some additional thresholding.
Neither of these methods guarantee that the pixels found are actually local maxima. Note where I said "larger or equal" above. Any plateau in your image (a region where all pixels have the same value) that is large enough will be marked by the algorithm, no matter if they are local maxima, local minima or somewhere in between.
For example:
import numpy as np
import matplotlib.pyplot as pp
import scipy.ndimage as ndimage
manos = np.sin(np.arange(100)/10)
manos = np.round(30*manos)/30     # Rounding to create plateaus
giorgos = ndimage.filters.maximum_filter(manos, footprint=np.ones([3]))
giorgos = (giorgos == manos)
pp.plot(manos);
pp.plot(giorgos);
pp.show()
Note how the filter identified three points near the local minimum of the sinusoid. The middle one of these is the actual local minimum, the other two are plateaus that are neither local maxima nor minima.
In contrast, the MATLAB function imregionalmax identifies all plateaus that are surrounded by pixels with a lower value. The algorithm required to do this is very different from the one above. It can be efficiently accomplished using a Union-Find algorithm, or less efficiently using a flood-fill-type algorithm. The main idea is to find a pixel that is not lower than any neighbor, then expand from it to its equal-valued neighbors until the whole plateau has been explored or until you find one of the pixels in the plateau with a higher-valued neighbor.
One implementation available from Python is in DIPlib (note: I'm an author):
import diplib as dip
nikos = dip.Maxima(manos)
pp.plot(manos);
pp.plot(nikos);
pp.show()
Another implementation is in SciKit-Image (Thanks to Juan for pointing this out):
nikos = skimage.morphology.local_maxima(manos)
                Note: we now have flood-fill based plateau finding in scikit-image as skimage.morphology.local_maxima
– Juan
                May 15, 2019 at 1:43
                Thanks!  I think local_maxima just works fine for me. @Juan, I got two future warnings in extrema.py lines 239 and 242. Any idea what this might be? thanks
– Hipparkhos
                May 15, 2019 at 10:42
                What version are you on? I think we fixed those warnings in github.com/scikit-image/scikit-image/pull/3238, specifically github.com/scikit-image/scikit-image/commit/… . If they come up in the latest version (0.14.2 or 0.15.0), please raise an issue on GitHub.
– Juan
                May 16, 2019 at 2:32
        
Thanks for contributing an answer to Stack Overflow!
Please be sure to answer the question. Provide details and share your research!
But avoid …
Asking for help, clarification, or responding to other answers.
Making statements based on opinion; back them up with references or personal experience.
To learn more, see our tips on writing great answers.