Numpy滤波器以平滑零区域

Question

我可以帮助您进行翻译。以下是需要翻译的内容：

我有一个由0及以上整数组成的2D numpy数组，其中数值代表区域标签。例如，

array([[9, 9, 9, 0, 0, 0, 0, 1, 1, 1],
       [9, 9, 9, 9, 0, 7, 1, 1, 1, 1],
       [9, 9, 9, 9, 0, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 0, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 0, 2, 2, 2, 1, 1],
       [4, 4, 4, 4, 0, 2, 2, 2, 1, 1],
       [4, 6, 6, 4, 0, 0, 0, 0, 0, 0],
       [4, 6, 6, 4, 0, 0, 0, 0, 0, 0],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5]])

我希望将指数等于0的值（即零区域）赋予其邻域中最常见的值。该操作本质上关闭了零区域。我尝试过多种膨胀、腐蚀、灰度闭合（grey-closing）和其他形态学操作（morphology operations）的变化，但是我无法完全消除零区域（而不是笨拙地混合其他区域）。一个不错的方法是定义一个仅在零上卷积的核，并在滤波区域中设置具有最常见标签的值。但我不确定如何实现这一点。

- BoltzmannBrain

---

你愿意使用Numba吗？ - John Zwinck

Answer 1

这里提出了一种向量化方法。步骤如下：

1.获取大小为2D的滑动窗口内核，导致4D数组。我们可以使用skimage的view_as_windows将其作为视图获取，从而避免为此创建任何额外的内存。

2.通过索引进入4D数组中心为零的窗口。这会强制进行复制。但是假设0的数量相对于输入数组中的总元素数量较小，则应该没问题。

3.对于每个选择的窗口，使用适当的偏移量来偏移每个窗口，以使用np.bincount执行计数。因此，使用bincount并获得排除零的最大计数。最大计数的argmax应该是我们所需的结果！

以下是覆盖这些步骤的实现：

from skimage.util import view_as_windows as viewW

def fill_zero_regions(a, kernel_size=3):
    hk = kernel_size//2 # half_kernel_size    

    a4D = viewW(a, (kernel_size,kernel_size))
    sliced_a = a[hk:-hk,hk:-hk]
    zeros_mask = sliced_a==0
    zero_neighs = a4D[zeros_mask].reshape(-1,kernel_size**2)
    n = len(zero_neighs) # num_zeros

    scale = zero_neighs.max()+1
    zno = zero_neighs + scale*np.arange(n)[:,None] # zero_neighs_offsetted

    count = np.bincount(zno.ravel(), minlength=n*scale).reshape(n,-1)
    modevals = count[:,1:].argmax(1)+1
    sliced_a[zeros_mask] = modevals
    return a

样例运行 -

In [23]: a
Out[23]: 
array([[9, 9, 9, 0, 0, 0, 0, 1, 1, 1],
       [9, 9, 9, 9, 0, 7, 1, 1, 1, 1],
       [9, 9, 9, 9, 0, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 0, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 0, 2, 2, 2, 1, 1],
       [4, 4, 4, 4, 0, 2, 2, 2, 1, 1],
       [4, 6, 6, 4, 0, 0, 0, 0, 0, 0],
       [4, 6, 6, 4, 0, 0, 0, 0, 0, 0],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5]])

In [24]: fill_zero_regions(a)
Out[24]: 
array([[9, 9, 9, 0, 0, 0, 0, 1, 1, 1],
       [9, 9, 9, 9, 9, 7, 1, 1, 1, 1],
       [9, 9, 9, 9, 2, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 2, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 2, 2, 2, 2, 1, 1],
       [4, 4, 4, 4, 2, 2, 2, 2, 1, 1],
       [4, 6, 6, 4, 4, 2, 2, 2, 1, 0],
       [4, 6, 6, 4, 4, 5, 5, 5, 5, 0],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5]])

正如所见，我们没有解决边界情况。如果需要，请使用零填充数组作为输入数组，类似于这样：np.pad(a, (k//2,k//2), 'constant')，其中k是内核大小（对于示例=3）。

- Divakar

---

非常好的回答，谢谢！+1 为解释内存使用。 - BoltzmannBrain

Answer 2

基于卷积思想的一种可能解决方案。

from scipy import stats
ar = #Your np array
blank = np.zeros(ar.shape)
#Size to search in for mode values
window_size = 3

for x,y in np.array(np.where(ar == 0)).T:
    window = ar[max(x-window_size,0):x+window_size,max(0,y-window_size):y+window_size]
    oneD = window.flatten()

    #fill blank array with modal value
    blank[x,y] = stats.mode(oneD[oneD != 0])[0]

#fill in the zeros
print ar + blank

我不确定在这里避免循环是否可能

- DJK

---

zip(np.where(ar == 0)[0],np.where(ar == 0)[1]) 可以简化为 np.array(np.where(ar == 0)).T。 - John Zwinck

Answer 3

这是一个使用Numba的工作解决方案，我尚未对其进行分析，但应该非常快速：

import numba
@numba.njit
def nn(arr):
    res = arr.copy()
    zeros = np.where(arr == 0)
    for n in range(len(zeros[0])):
        i = zeros[0][n]
        j = zeros[1][n]
        left = max(i-1, 0)
        right = min(i+2, arr.shape[1])
        top = max(j-1, 0)
        bottom = min(j+2, arr.shape[0])
        area = arr[left:right,top:bottom].ravel()
        counts = np.bincount(area[area != 0])
        res[i,j] = np.argmax(counts)
    return res

它产生：

array([[9, 9, 9, 9, 7, 1, 1, 1, 1, 1],
       [9, 9, 9, 9, 9, 7, 1, 1, 1, 1],
       [9, 9, 9, 9, 2, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 2, 2, 2, 1, 1, 1],
       [9, 9, 9, 8, 2, 2, 2, 2, 1, 1],
       [4, 4, 4, 4, 2, 2, 2, 2, 1, 1],
       [4, 6, 6, 4, 4, 2, 2, 2, 1, 1],
       [4, 6, 6, 4, 4, 5, 5, 5, 5, 5],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5],
       [4, 4, 4, 4, 5, 5, 5, 5, 5, 5]])

这里的核大小为3x3，由 i 和 j 减1并加2定义（加2是因为Python切片取一个超过末尾的元素，例如[0：3]给出3个元素）。边界条件由 min 和 max 处理。

bincount想法的来源：https://dev59.com/UW025IYBdhLWcg3wChSc#6252400

- John Zwinck