19

I have a numpy array like this:

x = np.array([[1,2,3],[4,5,6],[7,8,9]])

I need to create a function let's call it "neighbors" with the following input parameter:

  • x: a numpy 2d array
  • (i,j): the index of an element in a 2d array
  • d: neighborhood radius

As output I want to get the neighbors of the cell i,j with a given distance d. So if I run

neighbors(im, i, j, d=1) with i = 1 and j = 1 (element value = 5) 

I should get the indices of the following values: [1,2,3,4,6,7,8,9]. I hope I make it clear. Is there any library like scipy which deal with this?

I've done something working but it's a rough solution.

def pixel_neighbours(self, p):

    rows, cols = self.im.shape

    i, j = p[0], p[1]

    rmin = i - 1 if i - 1 >= 0 else 0
    rmax = i + 1 if i + 1 < rows else i

    cmin = j - 1 if j - 1 >= 0 else 0
    cmax = j + 1 if j + 1 < cols else j

    neighbours = []

    for x in xrange(rmin, rmax + 1):
        for y in xrange(cmin, cmax + 1):
            neighbours.append([x, y])
    neighbours.remove([p[0], p[1]])

    return neighbours

How can I improve this?

0

8 Answers 8

29

Have a look at scipy.ndimage.generic_filter.

As an example:

import numpy as np
import scipy.ndimage as ndimage

def test_func(values):
    print(values)
    return values.sum()


x = np.array([[1,2,3],[4,5,6],[7,8,9]])

footprint = np.array([[1,1,1],
                      [1,0,1],
                      [1,1,1]])

results = ndimage.generic_filter(x, test_func, footprint=footprint)

By default, it will "reflect" the values at the boundaries. You can control this with the mode keyword argument.

However, if you're wanting to do something like this, there's a good chance that you can express your problem as a convolution of some sort. If so, it will be much faster to break it down into convolutional steps and use more optimized functions (e.g. most of scipy.ndimage).

1
  • This method is slower than my custom own method using for loops
    – Angel
    Dec 31, 2021 at 12:21
10

EDIT: ah crap, my answer is just writing im[i-d:i+d+1, j-d:j+d+1].flatten() but written in a incomprehensible way :)


The good old sliding window trick may help here:

import numpy as np
from numpy.lib.stride_tricks import as_strided

def sliding_window(arr, window_size):
    """ Construct a sliding window view of the array"""
    arr = np.asarray(arr)
    window_size = int(window_size)
    if arr.ndim != 2:
        raise ValueError("need 2-D input")
    if not (window_size > 0):
        raise ValueError("need a positive window size")
    shape = (arr.shape[0] - window_size + 1,
             arr.shape[1] - window_size + 1,
             window_size, window_size)
    if shape[0] <= 0:
        shape = (1, shape[1], arr.shape[0], shape[3])
    if shape[1] <= 0:
        shape = (shape[0], 1, shape[2], arr.shape[1])
    strides = (arr.shape[1]*arr.itemsize, arr.itemsize,
               arr.shape[1]*arr.itemsize, arr.itemsize)
    return as_strided(arr, shape=shape, strides=strides)

def cell_neighbors(arr, i, j, d):
    """Return d-th neighbors of cell (i, j)"""
    w = sliding_window(arr, 2*d+1)

    ix = np.clip(i - d, 0, w.shape[0]-1)
    jx = np.clip(j - d, 0, w.shape[1]-1)

    i0 = max(0, i - d - ix)
    j0 = max(0, j - d - jx)
    i1 = w.shape[2] - max(0, d - i + ix)
    j1 = w.shape[3] - max(0, d - j + jx)

    return w[ix, jx][i0:i1,j0:j1].ravel()

x = np.arange(8*8).reshape(8, 8)
print x

for d in [1, 2]:
    for p in [(0,0), (0,1), (6,6), (8,8)]:
        print "-- d=%d, %r" % (d, p)
        print cell_neighbors(x, p[0], p[1], d=d)

Didn't do any timings here, but it's possible this version has reasonable performance.

For more info, search the net with phrases "rolling window numpy" or "sliding window numpy".

0
4

By using max and min, you handle pixels at the upper and lower bounds:

im[max(i-1,0):min(i+2,i_end), max(j-1,0):min(j+2,j_end)].flatten()
4

I don't know about any library functions for this, but you can easily write something like this yourself using the great slicing functionality of numpy:

import numpy as np
def neighbors(im, i, j, d=1):
    n = im[i-d:i+d+1, j-d:j+d+1].flatten()
    # remove the element (i,j)
    n = np.hstack((n[:len(n)//2], n[len(n)//2+1:] ))
    return n

Of course you should do some range checks to avoid out-of-bounds access.

2
  • 2
    What does b refer to here?
    – James
    Sep 13, 2017 at 10:26
  • b is not defined, as James said May 24, 2020 at 16:28
2

I agree with Joe Kingtons response, just an add to the footprints

import numpy as np
from scipy.ndimage import generate_binary_structure
from scipy.ndimage import iterate_structure
foot = np.array(generate_binary_structure(2, 1),dtype=int)

or for bigger/different footprints for ex.

np.array(iterate_structure(foot , 2),dtype=int)
0
from skimage.util import view_as_windows
window_shape = (3, 3, 3)
neighborhoods = view_as_windows(img[200:220], window_shape) 

Here is one D if you need:

>>> A = np.arange(10)
>>> A
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> window_shape = (3,)
>>> B = view_as_windows(A, window_shape)
>>> B.shape
(8, 3)
>>> B
array([[0, 1, 2],
       [1, 2, 3],
       [2, 3, 4],
       [3, 4, 5],
       [4, 5, 6],
       [5, 6, 7],
       [6, 7, 8],
       [7, 8, 9]]

reference :https://scikit-image.org/docs/stable/api/skimage.util.html#skimage.util.view_as_windows

-1

Possibly use a KDTree in SciPy ?

1
  • 2
    Just for whatever it's worth, for regularly-gridded data, a quadtree is not ideal. Indexing the grid that the data is on gives you a much faster lookup for neighbors. Jun 12, 2012 at 13:33
-1

We first init our matrix of interest using numpy

import numpy as np

x = np.array([[1,2,3],[4,5,6],[7,8,9]])

print(x)

[[1 2 3]
 [4 5 6]
 [7 8 9]]

Our neighbors is a function of distance for instance we might be interested in neighbors of distance 2 this tells us how should we pad our matrix x. We choose to pad with zeros but you can pad with whatever you like might be mean,mode,median of a row/column

d = 2

x_padded = np.pad(x,d,mode='constant')

print(x_padded)

[[0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0]
 [0 0 1 2 3 0 0]
 [0 0 4 5 6 0 0]
 [0 0 7 8 9 0 0]
 [0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0]]

We use x_padded matrix to get neighbors of any value in matrix x. Let (i,j) and (s,t) be indexes of x and x_padded respectively. Now we need to translate (i,j) to (s,t) to get neighbors of (i,j)

i,j = 2,1
s,t = 2*d+i+1, 2*d+j+1

window = x_padded[i:s, j:t]

print(window)

[[0 1 2 3 0]
 [0 4 5 6 0]
 [0 7 8 9 0]
 [0 0 0 0 0]
 [0 0 0 0 0]]

Please Note!!! the indexes (i,j) point to any value you wish to get its neighbors in matrix x

One might wish to iterate over each point in matrix x, get its neighbors and do computation using the neighbors for instance in Image Processing, the convolution with a kernel. One might do the following to get neighbors of each pixel in an image x

for i in range(x.shape[0]):
    for j in range(x.shape[1]):
        i,j = 2,1
        s,t = 2*d+i+1, 2*d+j+1
        window = x_padded[i:s, j:t]
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.