What's the difference between scipy.ndimage.filters.convolve and scipy.signal.convolve?

Question

As far as I have seen, these methods are both implemented as C functions in the respective DLLs, and it appears that the ndimage version is faster (neither implementation uses parallelized code, like calls to blas or MKL).

Also, when I tried to check that they return the same results by running the following code, the assertion of equality failed. I couldn't figure out from the documentation what exactly the functional differences between the two methods should be (the documentation isn't very clear about what 0 means relative to the location of the kernel's origin either; from examples, I deduced it's in the center, but I might be wrong).

from numpy import random, allclose
from scipy.ndimage.filters import convolve as convolveim
from scipy.signal import convolve as convolvesig

a = random.random((100, 100, 100))
b = random.random((10,10,10))

conv1 = convolveim(a,b, mode = 'constant')
conv2 = convolvesig(a,b, mode = 'same')

assert(allclose(conv1,conv2))

Thanks!

Answer 1

The two functions have different conventions for dealing with the boundary. To make your calls functionally the same, add the argument origin=-1 or origin=(-1,-1,-1) to the call to convolveim:

In [46]: a = random.random((100,100,100))

In [47]: b = random.random((10,10,10))

In [48]: c1 = convolveim(a, b, mode='constant', origin=-1)

In [49]: c2 = convolvesig(a, b, mode='same')

In [50]: allclose(c1,c2)
Out[50]: True

Only shift the origin when the dimensions of b are even. When they are odd, the functions agree when the default origin=0 is used:

In [88]: b = random.random((11,11,11))

In [89]: c1 = convolveim(a, b, mode='constant')

In [90]: c2 = convolvesig(a, b, mode='same')

In [91]: allclose(c1,c2)
Out[91]: True

Answer 2

There is a very important difference. The implementation in the image package seems to be the typical restricted version used in image processing to achieve the 'same' size of the image after the convolution. So it coincides with the 'same' option in the signal processing package if we use the mode='constant', as in the examples above. The signal processing package seems to implement the real strict definition of the convolution operator. Perhaps for this reason it is slower. Find enclosed some examples with completely different results.

In [13]: a=array([[1,2,1]])
In [14]: b=array([[1],[2],[1]])

In [17]: convolveim(a,b)
Out[17]: array([[4, 8, 4]])

In [18]: convolveim(b,a)
Out[18]: 
array([[4],
       [8],
       [4]])

In [19]: convolvesig(a,b)
Out[19]: 
array([[1, 2, 1],
       [2, 4, 2],
       [1, 2, 1]])

In [20]: convolvesig(b,a)
Out[20]: 
array([[1, 2, 1],
       [2, 4, 2],
       [1, 2, 1]])

Note that the signal processing package implementation is conmutative, as expected for a correct convolution. However, the implementation in the image package is not and it provides a solution with the 'same' dimensions as the first parameter.