CODEWITHSUNDEEP
pythonnumpymultidimensional-array
Aibar
Aibar

Reputation: 33

NumPy, how to effectively make elementwise operation involving other elements of ndarray(without loops)

Hi) I am trying to refactor my code with loops to NumPy operations to make code faster. Any clue how to do this? This code assigns a value to each element based on values of neighbor elements in 2D ndarray, and I couldn't find any answer for such specific staff.

This is the implementation of 6 neighbor method to find saddle points on photo described here https://documentcloud.adobe.com/link/track?uri=urn:aaid:scds:US:978c30d2-4888-491c-85c1-3949ea6166e9

It takes differences of current element with its neighbor elements. Then it counts changes of sign of those differences and if it >= 4 then it is saddle point.

Is it possible at all without loops?

Sorry if question is not clear or not in proper format - its my very first published question on StackOverflow

def findSaddlePoints6neibours(gray):
    gray = gray.astype(int)
    h = gray.shape[0]
    w = gray.shape[1]

    number = 0
    result = np.zeros((h, w))

    for y in range(1, h - 1):
        for x in range(1, w - 1):
            center = gray[y][x]
            neiboursDiff = []
            neiboursDiff.append(gray[y-1][x] - center)
            neiboursDiff.append(gray[y-1][x+1] - center)
            neiboursDiff.append(gray[y][x+1] - center)
            neiboursDiff.append(gray[y+1][x] - center)
            neiboursDiff.append(gray[y+1][x-1] - center)
            neiboursDiff.append(gray[y][x-1] - center)



            changes = 0
            for i in range(5):
                if (neiboursDiff[i] < 0 and neiboursDiff[i+1] > 0) or (neiboursDiff[i] > 0 and neiboursDiff[i+1] < 0):
                    changes += 1
            if (neiboursDiff[0] < 0 and neiboursDiff[5] > 0) or (neiboursDiff[0] > 0 and neiboursDiff[5] < 0):
                changes += 1
            if changes > 3:
                number += 1
                result[y][x] = 1

    return [result, number]

Upvotes: 3

Views: 122

Answers (1)

javidcf
javidcf

Reputation: 59701

Here is one vectorized solution:

import numpy as np

def findSaddlePoints6neibours_vec(gray):
    gray = np.asarray(gray, dtype=int)
    center = gray[1:-1, 1:-1]
    diffs = [
        gray[:-2, 1:-1],
        gray[:-2, 2:],
        gray[1:-1, 2:],
        gray[2:, 1:-1],
        gray[2:, :-2],
        gray[1:-1, :-2],
    ]
    diffs.append(diffs[0])
    diffs = np.stack(diffs)
    diffs -= center
    sign_changes = np.count_nonzero(diffs[:-1] * diffs[1:] < 0, axis=0)
    is_saddle = sign_changes > 3
    number = np.count_nonzero(is_saddle)
    result = np.pad(is_saddle, ((1, 1), (1, 1)), mode='constant').astype(int)
    return result, number

A quick test:

import numpy as np

# Make example input
np.random.seed(100)
gray = np.random.randint(-10, 10, size=(80, 100))

# The original function
result1, number1 = findSaddlePoints6neibours(gray)
# The vectorized function
result2, number2 = findSaddlePoints6neibours_vec(gray)
# Check results match
print(number1 == number2)
# True
print(np.all(result1 == result2))
# True

# Compare run times
%timeit findSaddlePoints6neibours(gray)
# 31.1 ms ± 682 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
%timeit findSaddlePoints6neibours_vec(gray)
# 247 µs ± 1.16 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

EDIT:

The disadvantage of the function above is that it takes more memory. If you can use Numba, you can compile the function and make it even faster if you use parallelization:

import numba as nb

@nb.njit(parallel=True)
def findSaddlePoints6neibours_nb(gray):
    gray = gray.astype(np.int32)
    h = gray.shape[0]
    w = gray.shape[1]
    number = 0
    result = np.zeros((h, w), dtype=np.int32)
    neiboursDiff = np.empty(7, dtype=np.int32)
    for y in nb.prange(1, h - 1):
        for x in np.prange(1, w - 1):
            neiboursDiff[0] = gray[y-1][x]
            neiboursDiff[1] = gray[y-1][x+1]
            neiboursDiff[2] = gray[y][x+1]
            neiboursDiff[3] = gray[y+1][x]
            neiboursDiff[4] = gray[y+1][x-1]
            neiboursDiff[5] = gray[y][x-1]
            neiboursDiff[6] = neiboursDiff[0]
            neiboursDiff -= gray[y, x]
            changes = np.sum(neiboursDiff[:-1] * neiboursDiff[1:] < 0)
            is_saddle = int(changes > 3)
            number += is_saddle
            result[y, x] = is_saddle
    return result, number

Continuing the small benchmark above:

%timeit findSaddlePoints6neibours_nb(gray)
# 114 µs ± 496 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

Upvotes: 2

Related Questions