Convolve2d 只需使用 Numpy

Posted 2023-02-23

【中文标题】Convolve2d 只需使用 Numpy

【英文标题】：Convolve2d just by using Numpy

【发布时间】：2017-08-22 12:44:09

【问题描述】：

我正在研究使用 Numpy 进行图像处理，但遇到了卷积过滤的问题。

我想对灰度图像进行卷积。 （将二维数组与较小的二维数组卷积）

有人想改进我的方法吗？

我知道scipy 支持 convolve2d，但我只想使用 Numpy 制作 convolve2d。

我做了什么

首先，我为子矩阵制作了一个二维数组。

a = np.arange(25).reshape(5,5) # original matrix

submatrices = np.array([
     [a[:-2,:-2], a[:-2,1:-1], a[:-2,2:]],
     [a[1:-1,:-2], a[1:-1,1:-1], a[1:-1,2:]],
     [a[2:,:-2], a[2:,1:-1], a[2:,2:]]])

子矩阵看起来很复杂，但我正在做的事情如下图所示。

接下来，我将每个子矩阵与一个过滤器相乘。

conv_filter = np.array([[0,-1,0],[-1,4,-1],[0,-1,0]])
multiplied_subs = np.einsum('ij,ijkl->ijkl',conv_filter,submatrices)

并将它们相加。

np.sum(np.sum(multiplied_subs, axis = -3), axis = -3)
#array([[ 6,  7,  8],
#       [11, 12, 13],
#       [16, 17, 18]])

因此这个过程可以称为我的 convolve2d。

def my_convolve2d(a, conv_filter):
    submatrices = np.array([
         [a[:-2,:-2], a[:-2,1:-1], a[:-2,2:]],
         [a[1:-1,:-2], a[1:-1,1:-1], a[1:-1,2:]],
         [a[2:,:-2], a[2:,1:-1], a[2:,2:]]])
    multiplied_subs = np.einsum('ij,ijkl->ijkl',conv_filter,submatrices)
    return np.sum(np.sum(multiplied_subs, axis = -3), axis = -3)

但是，我发现这个 my_convolve2d 很麻烦，原因有 3 个。

子矩阵的生成过于笨拙，难以阅读，只能在过滤器为3*3时使用 各种子矩阵的大小似乎太大了，因为它比原始矩阵大了大约 9 倍。 求和似乎有点不直观。简单地说，丑。

感谢您阅读本文。

某种更新。我为自己写了一个conv3d。我会将其保留为公共领域。

def convolve3d(img, kernel):
    # calc the size of the array of submatracies
    sub_shape = tuple(np.subtract(img.shape, kernel.shape) + 1)

    # alias for the function
    strd = np.lib.stride_tricks.as_strided

    # make an array of submatracies
    submatrices = strd(img,kernel.shape + sub_shape,img.strides * 2)

    # sum the submatraces and kernel
    convolved_matrix = np.einsum('hij,hijklm->klm', kernel, submatrices)

    return convolved_matrix

【问题讨论】：

感谢您提供矩阵图 :) 如果我理解正确，您需要有关如何使您的解决方案更优雅的提示吗？

很高兴它有帮助！是的。如果您能提供一些技巧来克服最后几行中写的 3 个问题，我将不胜感激。

我应该补充一点，这 3 个点是按优先顺序排列的。第一个对我来说很重要，最后一个似乎有点微不足道。如果还有其他问题和改进，我也会很高兴。

第二张图（等号后）不是错了吗？不应该将每个子矩阵与过滤器相乘（按元素），然后将每个结果子矩阵的元素相加吗？

@AndyK 他们会产生同样的结果。

【参考方案1】：

您可以使用as_strided 生成子数组：

import numpy as np

a = np.array([[ 0,  1,  2,  3,  4],
       [ 5,  6,  7,  8,  9],
       [10, 11, 12, 13, 14],
       [15, 16, 17, 18, 19],
       [20, 21, 22, 23, 24]])

sub_shape = (3,3)
view_shape = tuple(np.subtract(a.shape, sub_shape) + 1) + sub_shape
strides = a.strides + a.strides

sub_matrices = np.lib.stride_tricks.as_strided(a,view_shape,strides)

要摆脱第二个“丑陋”的总和，请更改 einsum，使输出数组只有 j 和 k。这意味着您的第二次总结。

conv_filter = np.array([[0,-1,0],[-1,5,-1],[0,-1,0]])
m = np.einsum('ij,ijkl->kl',conv_filter,sub_matrices)

# [[ 6  7  8]
#  [11 12 13]
#  [16 17 18]]

【讨论】：

如果 a_s 是跨步数组并且 filter 是您的类似拉普拉斯算子的过滤器，那么如果您的答案确实是 array( [[ 6, 7, 8], [11, 12, 13], [16, 17, 18]])

感谢您的提示。我现在正在自己尝试这个。也许微不足道，但我认为名称过滤器不合适，因为它是python的内置函数。

可以直接在爱因斯坦求和中进行求和。查看答案

这个问题，写sub_shape = conv_filter.shape不是更好吗？

看来不一定是这样。 conv_filter.shape 应该是 a.shape - sub_shape + 1。在您的示例a.shape = [5, 5] 和sub_shape = (3, 3) 中就是这样，所以conv_filter.shape == sub_shape 是正确的，但没有严格执行。如果你想要它被强制执行，那么你应该设置sub_shape = a.shape - conv_filter.shape + 1

【参考方案2】：

使用as_strided 和@Crispin 的einsum 技巧从上方清理。将过滤器尺寸强制为扩展形状。如果索引兼容，甚至应该允许非方形输入。

def conv2d(a, f):
    s = f.shape + tuple(np.subtract(a.shape, f.shape) + 1)
    strd = numpy.lib.stride_tricks.as_strided
    subM = strd(a, shape = s, strides = a.strides * 2)
    return np.einsum('ij,ijkl->kl', f, subM)

【讨论】：

进一步简化...请参阅下面的评论... np.sum(a_s * filter, axis=(2,3)) 如果您的答案确实是 array([[ 6, 7, 8 ], [11, 12, 13], [16, 17, 18]]) ... 其中 a_s 是跨步数组，filter 是 3x3 过滤器

不确定为什么会这样@NaN，因为它肯定没有按照问题的要求做 - 但它确实做到了，即使对于任意的a 矩阵

至少在 numpy v12 中，a.shape 和 f.shape 是一个元组，所以我认为 s 应该是 tuple(np.subtract(a.shape, f.shape)+1)。

真的！仍然习惯于shape 给出ndarray。解决了这个问题。

嗨@DanielF，是否有适用于RGB的概括？对 einsum 表示法不太熟悉，所以知道如何概括会很棒。

【参考方案3】：

您还可以使用 fft（执行卷积的更快方法之一）

from numpy.fft import fft2, ifft2
import numpy as np

def fft_convolve2d(x,y):
    """ 2D convolution, using FFT"""
    fr = fft2(x)
    fr2 = fft2(np.flipud(np.fliplr(y)))
    m,n = fr.shape
    cc = np.real(ifft2(fr*fr2))
    cc = np.roll(cc, -m/2+1,axis=0)
    cc = np.roll(cc, -n/2+1,axis=1)
    return cc

https://gist.github.com/thearn/5424195 您必须将滤镜填充为与图像相同的大小（将其放在 zeros_like 垫子的中间。）

干杯， 丹

【讨论】：

如果内核只有 9 个元素，使用 FFT 进行卷积肯定效率不高。这仅对大内核有用。

是的。它取决于内核和图像大小，但 fft 表现出色的阈值非常低。它还取决于您的计算机架构，但对于大多数通用用例来说，fft 是高性能库的首选。

下面的文章比较了各种图像卷积方法的性能，由此我们可以得出结论，Numpy的FFT方法是执行卷积最快的方法。 laurentperrinet.github.io/sciblog/posts/…

【参考方案4】：

https://laurentperrinet.github.io/sciblog/posts/2017-09-20-the-fastest-2d-convolution-in-the-world.html

在此处查看所有卷积方法及其各自的性能。 另外，我发现下面的代码 sn-p 更简单。

from numpy.fft  import fft2, ifft2
def np_fftconvolve(A, B):
    return np.real(ifft2(fft2(A)*fft2(B, s=A.shape)))

【讨论】：

此代码不正确。它将改变输出图像。内核越大，偏移越大。

你能详细说明它将如何改变输出图像吗？

***.com/questions/54877892/…