2
\$\begingroup\$

I'm trying to find the local maxima in a 3D NumPy array, but I can't seem to find an easy way to do that using NumPy, SciPy, or anything else.

For now I implemented it using scipy.signal.argrelextrema.
But it takes very long to process large arrays, and only works on separated axes.

import numpy as np
from scipy.signal import argrelextrema


def local_maxima_3D(data, order=1):
    """Detects local maxima in a 3D array

    Parameters
    ---------
    data : 3d ndarray
    order : int
        How many points on each side to use for the comparison

    Returns
    -------
    coordinates : ndarray
        coordinates of the local maxima
    values : ndarray
        values of the local maxima
    """
    # Coordinates of local maxima along each axis
    peaks0 = np.array(argrelextrema(data, np.greater, axis=0, order=order))
    peaks1 = np.array(argrelextrema(data, np.greater, axis=1, order=order))
    peaks2 = np.array(argrelextrema(data, np.greater, axis=2, order=order))
    
    # Stack all coordinates
    stacked = np.vstack((peaks0.transpose(), peaks1.transpose(),
                         peaks2.transpose()))

    # We keep coordinates that appear three times (once for each axis)
    elements, counts = np.unique(stacked, axis=0, return_counts=True)
    coords = elements[np.where(counts == 3)[0]]

    # Compute values at filtered coordinates
    values = data[coords[:, 0], coords[:, 1], coords[:, 2]]

    return coords, values

I know this solution is far from optimal and only works with order=1.
Is there any better way to find local maxima in a 3D array in Python?

\$\endgroup\$
3
  • 3
    \$\begingroup\$ Do you know anything about the data? If they're determined by a formula, then you should consider getting an analytical answer. Otherwise, the data being continuous, derivable and/or single-maximum would allow for your use of an optimizer, especially if you can provide a gradient function. \$\endgroup\$ Commented Apr 1, 2019 at 13:54
  • \$\begingroup\$ As opposed to unique(stacked(...))[where(counts==3)], have you tried set intersection, like set(peaks0) & set(peaks1) & set(peaks2). Of course, the peaks# coordinates would need to be converted into something hashable (tuples) in order to create the sets. I don’t know if this would be faster, so this is an experiment not an answer ... \$\endgroup\$ Commented Apr 1, 2019 at 14:28
  • \$\begingroup\$ It would be better to put down, in the doc comment, what exactly constitutes a local maximum. \$\endgroup\$ Commented Jun 13 at 14:25

0

Reset to default

You must log in to answer this question.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.