I have a 3d numpy array and want to generate a secondary array consisting of the minimum of each value and the values in the 10 rows directly above and 10 rows directly below (i.e each entry is the minimum value from 21 values) for each 2d array.
I've been trying to use 'numpy.clip' to deal with the edges of the array - here the range of values which the minimum is taken from should simply reduce to 10 at the values on the top/bottom of the array. I think something like 'scipy.signal.argrelmin' seems to be what I'm after.
Here's my code so far, definitely not the best way to go about it:
import numpy as np
array_3d = np.random.random_integers(50, 80, (3, 50, 18))
minimums = np.zeros(array_3d.shape)
for array_2d_index in range(len(array_3d)):
for row_index in range(len(array_3d[array_2d_index])):
for col_index in range(len(array_3d[array_2d_index][row_index])):
minimums[array_2d_index][row_index][col_index] = min(array_3d[array_2d_index][np.clip(row_index-10, 0, 49):np.clip(row_index+10, 0, 49)][col_index])
The main issue I think is that this is taking the minimum from the columns either side of each entry instead of the rows, which has been giving index errors.
Would appreciate any advice, thanks.
Approach #1
Here's one approach with np.lib.stride_tricks.as_strided -
def strided_3D_axis1(array_3d, L):
s0,s1,s2 = array_3d.strides
strided = np.lib.stride_tricks.as_strided
m,n,r = array_3d.shape
nL = n-L+1
return strided(array_3d, (m,nL,L,r),(s0,s1,s1,s2))
out = strided_3D_axis1(array_3d, L=21).min(axis=-2)
Sample run -
1) Input :
In [179]: array_3d
Out[179]:
array([[[73, 65, 51, 76, 59],
[74, 57, 75, 53, 70],
[60, 74, 52, 54, 60],
[54, 52, 62, 75, 50],
[68, 56, 68, 63, 77]],
[[62, 70, 60, 79, 74],
[70, 68, 50, 74, 57],
[63, 57, 69, 65, 54],
[63, 63, 68, 58, 60],
[70, 66, 65, 78, 78]]])
2) Strided view :
In [180]: strided_3D_axis1(array_3d, L=3)
Out[180]:
array([[[[73, 65, 51, 76, 59],
[74, 57, 75, 53, 70],
[60, 74, 52, 54, 60]],
[[74, 57, 75, 53, 70],
[60, 74, 52, 54, 60],
[54, 52, 62, 75, 50]],
[[60, 74, 52, 54, 60],
[54, 52, 62, 75, 50],
[68, 56, 68, 63, 77]]],
[[[62, 70, 60, 79, 74],
[70, 68, 50, 74, 57],
[63, 57, 69, 65, 54]],
[[70, 68, 50, 74, 57],
[63, 57, 69, 65, 54],
[63, 63, 68, 58, 60]],
[[63, 57, 69, 65, 54],
[63, 63, 68, 58, 60],
[70, 66, 65, 78, 78]]]])
3) Strided view based min :
In [181]: strided_3D_axis1(array_3d, L=3).min(axis=-2)
Out[181]:
array([[[60, 57, 51, 53, 59],
[54, 52, 52, 53, 50],
[54, 52, 52, 54, 50]],
[[62, 57, 50, 65, 54],
[63, 57, 50, 58, 54],
[63, 57, 65, 58, 54]]])
Approach #2
Here's another with broadcasting upon creating all sliding indices along the second axis -
array_3d[:,np.arange(array_3d.shape[1]-L+1)[:,None] + range(L)].min(-2)
Approach #3
Here's another using Scipy's 1D minimum filter -
from scipy.ndimage.filters import minimum_filter1d as minf
L = 21
hL = (L-1)//2
out = minf(array_3d,L,axis=1)[:,hL:-hL]
Runtime test -
In [231]: array_3d = np.random.randint(50, 80, (3, 50, 18))
In [232]: %timeit strided_3D_axis1(array_3d, L=21).min(axis=-2)
10000 loops, best of 3: 54.2 µs per loop
In [233]: %timeit array_3d[:,np.arange(array_3d.shape[1]-L+1)[:,None] + range(L)].min(-2)
10000 loops, best of 3: 81.3 µs per loop
In [234]: L = 21
...: hL = (L-1)//2
...:
In [235]: %timeit minf(array_3d,L,axis=1)[:,hL:-hL]
10000 loops, best of 3: 32 µs per loop