1

I have a numpy array a:

a = np.array([[0,4,3,9,9,9,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,4,3,5,9,9,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,8,3,9,2,6,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,9,2,6,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,4,2,9,2,6,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,2,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,9,2,6,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,2,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,4,2,2,6,8,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,4,4,3,4,4,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10],
          [4,4,4,3,4,4,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10]])

I'm looking for the index of the smallest neighbours of cell(s) where value= 2 (numpy.where(a==value)) but bigger than value. I also need the index of the corresponding cell(s) = value for which we found the smallest neighbour.

Result in this case (value = 2) should be:

  • index of neighbouring cells: [0,2] and [4,3] The corresponding and

  • index of corresponding cell(s): [1,2] and [3,3]

Please apologies if the question is not very clear.

This is what I have so far: 

import numpy as np
value = 2
neighbors = np.zeros(4, dtype = np.float)
fdx = np.flatnonzero(a== value)
locations = fdx // a.shape[1], fdx % a.shape[1]
maximums = []
for item in zip(*locations):
    i, j = item[0], item[1]
    neighbors[0], neighbors[1], neighbors[2],neighbors[3] = [a[i-1,j], a[i+1,j], a[i,j-1], a[i,j+1]]
    maximums.append(min(neighbors[neighbors> value]))
print np.where(a==min(maximums))

prints: (array([0, 4]), array([2, 3]))

Which is very slow and also still don't know how to find the index of corresponding cells. Any solution which is totally different to my solution would be also accepted.

Behzad Jamali
  • 884
  • 2
  • 10
  • 23
  • Does `value` ever occur on the edges? If so your code will fail – Daniel F Nov 28 '17 at 07:48
  • @DanielF Yes it does! I didn't ask for considering edges but it is indeed important. Your solution comes at expense of lower speed though (in comparison to the answer provided by Paul Panzer and Antoine Zambelli (I just looped each solution 10000 times to compare their speed which might not be a good way to compare them). However, I guess speed comparison result might be different if tested for much larger arrays? – Behzad Jamali Nov 28 '17 at 07:55
  • Yeah, `binary_dilation` will probably be faster as it doesn't require allocating extra memory for the padded array. That's why I asked about the edges - if I didn't have to `pad`, using windowed views would probably be faster. – Daniel F Nov 28 '17 at 08:00
  • I wish I could accept multiple answers in SO! – Behzad Jamali Nov 28 '17 at 08:06
  • @DanielF I actually checked all the three answers for the edges and surprisingly they all work without a problem! Am I missing sth? tried `value = 10` and `value = 0` – Behzad Jamali Nov 28 '17 at 08:53
  • No, you're not missing anything. They all work on the edges. My statement was that my answer could probably be the fastest if you didn't *need* to account for the edges (i.e. I didn't need to `pad`) – Daniel F Nov 28 '17 at 10:14

3 Answers3

2

You can find neighbors using scipy.ndimage.morphology.binary_dilation

import numpy as np
from scipy.ndimage import morphology
a = np.array([[0,4,3,9,9,9],
              [4,4,2,2,2,9],
              [4,2,2,9,2,6],
              [4,2,2,2,6,8],
              [4,4,4,3,4,4]])
k = 2

# make mask
eq = a==k
# find greater neighbors (as mask)
nbs = morphology.binary_dilation(eq) & (a>k)
# translate to index
minidx = np.argwhere(nbs & (a == np.min(a[nbs])))

# now find neighbors' neighbors
# pad original mask
m,n = a.shape
eqp = np.zeros((m+2, n+2), bool)
eqp[1:-1,1:-1] = eq
# generate offset vectors for the four major directions (up, down, left, right)
# corrected for padding
offsp = np.array([(0,1),(2,1),(1,0),(1,2)])
# without padding correction
offs = offsp - 1#
# for each minimum, find all (1-4) reference neighbors
refidx = [i + offs[eqp[tuple((i+offsp).T)]] for i in minidx]

print(minidx)
print(refidx)

# [[0 2]
#  [4 3]]
# [array([[1, 2]]), array([[3, 3]])]
Paul Panzer
  • 51,835
  • 3
  • 54
  • 99
1

Not sure how much faster (if at all) this is, but it finds the corresponding cells.

import numpy as np

a = np.array([[0,4,3,9,9,9],
          [4,4,2,2,2,9],
          [4,2,2,9,2,6],
          [4,2,2,2,6,8],
          [4,4,4,3,4,4]])

value = 2
idx = np.where(a == value)
idx = zip(idx[0],idx[1])

running_min = np.inf
corresponding_idx = []#corresponding cells
nb_min_list = []#location of neighbors
nb_min_idx = []
for i in idx:
    nb_idx = [(i[0]+1,i[1]),(i[0]-1,i[1]),(i[0],i[1]+1),(i[0],i[1]-1)]#note no check for out of bounds.
    nb_idx = [nb for nb in nb_idx if nb[0] >= 0 and nb[0] < a.shape[0] and nb[1] >= 0 and nb[1] < a.shape[1]]#test for edges

    try:
        nb_min = min([a[nb] for nb in nb_idx if a[nb] > value])
        corresponding_idx.append(i)
        nb_min_list.append(nb_min)
        nb_min_idx.append([nb for nb in nb_idx if a[nb] == nb_min])
    except:
        pass

nb_min_loc = np.where(nb_min_list == min(nb_min_list))[0]
corresponding_cells = []
min_nbs = []
for nb in nb_min_loc:
    corresponding_cells.append(corresponding_idx[nb])
    min_nbs.append(nb_min_idx[nb])

print(corresponding_cells)#[(1, 2), (3, 3)]
print(min_nbs)#[[(0, 2)], [(4, 3)]]
Antoine Zambelli
  • 724
  • 7
  • 19
  • this solution was faster in a 10000 run loop. But it had an error while testing for another array: `nb_min = min([a[nb] for nb in nb_idx if a[nb] > value])` `ValueError: min() arg is an empty sequence` – Behzad Jamali Nov 28 '17 at 08:00
  • Fixed in the edit! It's a behavior of min()/max() for empty lists - it must have occurred in a situation where a cell with value 2 had all neighbors =2 (so the list came out empty. – Antoine Zambelli Nov 28 '17 at 08:19
  • Thanks! It works now! A bit slower than @Paul Penzer's answer though. – Behzad Jamali Nov 28 '17 at 08:35
  • The speed is very dependent on the `value` we choose. When the number of matching items with `value` is small, this solution is the fastest, otherwise it becomes very slow. – Behzad Jamali Nov 28 '17 at 09:04
  • Glad to hear it works. Interesting that it is so sensitive to 'sparsity' of the data (relative to `value`), I wouldn't have guessed. I'd go with Paul's answer then! – Antoine Zambelli Nov 28 '17 at 09:09
1

Using window_nd from a previous answer

def min_search(a, val = 2):
        a_view = window_nd(np.pad(a, ((1, 1),(1, 1)), 'constant',
                           constant_values = np.inf), 3)
        min_val = np.where(a_view[a == val] <= val, np.inf, a_view[a == val]).min()
        neig_mask = np.array([[0, 1, 0], [1, 0, 1], [0, 1, 0]], dtype = bool)
        rev_mask = np.logical_and(np.any(np.logical_and(a_view == min_val,
                                         neig_mask), axis = (-1, -2)), a == val)
        min_mask = np.logical_and(np.any(np.logical_and(a_view == val,
                                         neig_mask), axis = (-1, -2)), a == min_val)

        return np.nonzero(min_mask), np.nonzero(rev_mask)

What this does:

  • Creates a sliding window over a padded with val (so the windowed shape is (*a.shape, 3, 3))
  • finds the minimum value greater than val among all windows with val at the center (padding allows inlcusion of edges) an assigns it to min_val
  • neig_mask restricts neigbors to cardinal directions
  • finds windows with min_val in neig_mask positions and val in the center, assigns to rev_mask
  • finds windows with val in neig_mask positions and min_val in the center, assigns to min_mask

  • returns np.nonzero of min_mask and rev_mask, which are tuples that can be used as an index on a

    min_search(a)
Out:
((array([0, 4], dtype=int32), array([2, 3], dtype=int32)),
(array([1, 3], dtype=int32), array([2, 3], dtype=int32)))
Daniel F
  • 13,620
  • 2
  • 29
  • 55
  • Just one problem. I changed the example array (edited in the question) and realised that your solution returns all the cells that are equal to 3. I need only those cells that are in the neighbourhood of the selected cells. Please check your solution for the new array. – Behzad Jamali Nov 28 '17 at 08:15
  • yeah, turns out padding with `val` was a bad idea – Daniel F Nov 28 '17 at 08:35