9

I am trying to find a vectorized approach of finding the first position in an array where the values did not get higher than the maximum of n previous numbers. I thought about using the find_peaks method of scipy.signal to find a local maximum. I think it does exactly that if you define the distance to let's say 10 n is 10. But unfortunately, the condition for the distance has to be fulfilled in both directions - previous and upcoming numbers. Is there any other method or approach to finding such a thing?

Example:

arr1 = np.array([1.        , 0.73381293, 0.75649351, 0.77693474, 0.77884614,
       0.81055903, 0.81402439, 0.78798586, 0.78839588, 0.82967961,
       0.8448    , 0.83276451, 0.82539684, 0.81762916, 0.82722515,
       0.82101804, 0.82871127, 0.82825041, 0.82086957, 0.8347826 ,
       0.82666665, 0.82352942, 0.81270903, 0.81191224, 0.83180428,
       0.84975767, 0.84044236, 0.85057473, 0.8394649 , 0.80000001,
       0.83870965, 0.83962262, 0.85039371, 0.83359748, 0.84019768,
       0.83281732, 0.83660132])

from scipy.signal import find_peaks
peaks, _ = find_peaks(arr1, distance=10)

In this case, it finds positions 10 and 27. But also position 0 has 10 following elements which are not higher. How can I find those?

CDJB
  • 14,043
  • 5
  • 29
  • 55
Varlor
  • 1,421
  • 3
  • 22
  • 46

3 Answers3

6

Unfortunately, find_peaks() works by comparing neighboring values - so will not identify peaks that occur at the beginning or end of the array. One workaround is to use np.concatenate() to insert the minimum value of the array at the beginning and end, and then subtract 1 from the peaks variable:

>>> import numpy as np
>>> peaks, _ = find_peaks(np.concatenate(([min(arr1)],arr1,[min(arr1)])), distance=10)
>>> peaks-1
array([ 0, 10, 27], dtype=int64)
CDJB
  • 14,043
  • 5
  • 29
  • 55
2
def rolling_window(a, window):
    shape = a.shape[:-1] + (a.shape[-1] - window + 1, window)
    strides = a.strides + (a.strides[-1],)
    return np.lib.stride_tricks.as_strided(a, shape=shape, strides=strides)

def get_peaks(arr, window):
    maxss = np.argmax(rolling_window(arr1, window), axis=1)
    return np.where(maxss == 0)[0]
>>> arr1 = np.array([1.        , 0.73381293, 0.75649351, 0.77693474, 0.77884614,
       0.81055903, 0.81402439, 0.78798586, 0.78839588, 0.82967961,
       0.8448    , 0.83276451, 0.82539684, 0.81762916, 0.82722515,
       0.82101804, 0.82871127, 0.82825041, 0.82086957, 0.8347826 ,
       0.82666665, 0.82352942, 0.81270903, 0.81191224, 0.83180428,
       0.84975767, 0.84044236, 0.85057473, 0.8394649 , 0.80000001,
       0.83870965, 0.83962262, 0.85039371, 0.83359748, 0.84019768,
       0.83281732, 0.83660132])
>>> get_peaks(arr1, 10)
array([ 0, 10, 27])

Credit for rolling window function : Rolling window for 1D arrays in Numpy?

Frederik Bode
  • 2,632
  • 1
  • 10
  • 17
1

We can use 1D sliding-windowed max filter from SciPy. Also, it seems you are comparing against n previous elements. Since, the first element won't have any previous element, we need to make it ignore.

Hence, we would have the implementation, like so -

from scipy.ndimage.filters import maximum_filter1d

def peaks_previousN(a, n):
    W = (n-1)//2
    return np.flatnonzero(a[1:]>=maximum_filter1d(a, n, origin=W)[:-1])+1

Sample run with given sample array -

In [107]: peaks_previousN(arr1, n=10)
Out[107]: array([25, 27])
Divakar
  • 218,885
  • 19
  • 262
  • 358
  • You are right. Maybe i got confused when looking for a right way to do it. Im looking for the first position where the maximum did not changed 10 elements long. So it would be index 10. Because the maximum of 1 did not changed for 10 elements. So i think find peaks of scipy signal is exactly what im looking for. – Varlor Feb 06 '20 at 14:05
  • @Varlor So, are you saying compare the current against the last 10 element including or excluding the current element? – Divakar Feb 06 '20 at 14:08
  • The current against the last 10, excluding the current element. And i just need the first occurence – Varlor Feb 06 '20 at 14:11
  • @Varlor Then, I think the final expected output for `n = 10`, would be `[25,27]`. `10` won't be there because `arr1[0:10].max()` is `1.0` and `a[10]` is `0.8448`. Again, `0` won't be there because there are no previous elements to it. Am I right? – Divakar Feb 06 '20 at 14:18
  • 0 wouldnt be there exactly, because we are looking for the index where the overall maximum did not change n elements long. But it would be 10, because a[10] is 0.8448 and not higher than 1. So the maximum is still 1.0 and fullfiles the condition. Sorry if i was not clear enough prevously – Varlor Feb 06 '20 at 14:24
  • @Varlor Thought you were looking for peak. If `a[10]`, which is `0.8448 ` is not higher or equal to the max of previous 10 elements, how can it be called a peak? – Divakar Feb 06 '20 at 14:27
  • you are absoulutley right. In my way to search a vectorized solution of my problem solution i just thought to use somehow the peak function. It was a failure of me – Varlor Feb 06 '20 at 14:33