提高代码效率：滑动窗口上的标准差

Question

我正在尝试改进函数，该函数为图像的每个像素计算位于像素附近的像素的标准偏差。我的函数使用两个嵌入的循环来遍历矩阵，这是我的程序的瓶颈。我想可能有一种方法可以通过去掉numpy的循环来改进它，但是我不知道如何继续。 欢迎任何建议！

问候

def sliding_std_dev(image_original,radius=5) :
    height, width = image_original.shape
    result = np.zeros_like(image_original) # initialize the output matrix
    hgt = range(radius,height-radius)
    wdt = range(radius,width-radius)
    for i in hgt:
        for j in wdt:
            result[i,j] = np.std(image_original[i-radius:i+radius,j-radius:j+radius])
    return result

Answer 1

在图像处理中，最常用的方法是使用求和面积表，这是1984年this paper中引入的一种思想。其思想是，当你通过在一个窗口上加上来计算一个数量，并将窗口（例如，向右移动一个像素）时，你不需要添加新窗口中的所有项，你只需要从总数中减去最左边的列，然后添加新的最右边的列。因此，如果从数组中在两个维度上创建一个累加和数组，则可以在一个窗口上通过两个和和和一个减法得到和。如果为数组及其平方保留面积求和表，那么很容易从这两个表中得到方差。下面是一个实现：

def windowed_sum(a, win):
    table = np.cumsum(np.cumsum(a, axis=0), axis=1)
    win_sum = np.empty(tuple(np.subtract(a.shape, win-1)))
    win_sum[0,0] = table[win-1, win-1]
    win_sum[0, 1:] = table[win-1, win:] - table[win-1, :-win]
    win_sum[1:, 0] = table[win:, win-1] - table[:-win, win-1]
    win_sum[1:, 1:] = (table[win:, win:] + table[:-win, :-win] -
                       table[win:, :-win] - table[:-win, win:])
    return win_sum

def windowed_var(a, win):
    win_a = windowed_sum(a, win)
    win_a2 = windowed_sum(a*a, win)
    return (win_a2 - win_a * win_a / win/ win) / win / win

要确保这一点：

>>> a = np.arange(25).reshape(5,5)
>>> windowed_var(a, 3)
array([[ 17.33333333,  17.33333333,  17.33333333],
       [ 17.33333333,  17.33333333,  17.33333333],
       [ 17.33333333,  17.33333333,  17.33333333]])
>>> np.var(a[:3, :3])
17.333333333333332
>>> np.var(a[-3:, -3:])
17.333333333333332

这应该比基于卷积的方法快一些。