Taking an outer subtraction between a list of tuples using Numpy

Question

NOTE: I am not looking for syntax on how to use np.subtract.outer, instead I am trying to solve a very specific issue that I am facing during its application.

I have a list of tuples -

a = [(0,0), (1,0), (1,1), (2,0), (2,2)]
a = np.array(a)

I am trying to implement a full vectorized subtraction between each element of this list (tuples). I can do it with iterating using itertools as shown below expecting a 5x5x2 as -

#My attempt with iteration (producing expected results)

np.array([np.array(i[0])-np.array(i[1]) for i in itertools.product(a, a)]).reshape(5,5,2)

#RESULT (5x5x2) 
array([[[ 0,  0],
        [-1,  0],
        [-1, -1],
        [-2,  0],
        [-2, -2]],

       [[ 1,  0],
        [ 0,  0],
        [ 0, -1],
        [-1,  0],
        [-1, -2]],

       [[ 1,  1],
        [ 0,  1],
        [ 0,  0],
        [-1,  1],
        [-1, -1]],

       [[ 2,  0],
        [ 1,  0],
        [ 1, -1],
        [ 0,  0],
        [ 0, -2]],

       [[ 2,  2],
        [ 1,  2],
        [ 1,  1],
        [ 0,  2],
        [ 0,  0]]])

However, when using np.subtract.outer(). I end up with an array of shape 5x2x5x2 as below -

#My attempt with vectorization (producing unexpected results)

np.subtract.outer(a,a) #5x2x5x2

array([[[[ 0,  0],  #This matrix is same as ....
         [-1,  0],
         [-1, -1],
         [-2,  0],
         [-2, -2]],

        [[ 0,  0],  #This one...
         [-1,  0],
         [-1, -1],
         [-2,  0],
         [-2, -2]]],


       [[[ 1,  1],   #However, this ...
         [ 0,  1],
         [ 0,  0],
         [-1,  1],
         [-1, -1]],

        [[ 0,  0],  #And this, dont repeat! ..
         [-1,  0],
         [-1, -1],
         [-2,  0],
         [-2, -2]]],


       [[[ 1,  1],   #This ...
         [ 0,  1],
         [ 0,  0],
         [-1,  1],
         [-1, -1]],

        [[ 1,  1],   #And this do..
         [ 0,  1],
         [ 0,  0],
         [-1,  1],
         [-1, -1]]],


       [[[ 2,  2],
         [ 1,  2],
         [ 1,  1],
         [ 0,  2],
         [ 0,  0]],

        [[ 0,  0],
         [-1,  0],
         [-1, -1],
         [-2,  0],
         [-2, -2]]],


       [[[ 2,  2],
         [ 1,  2],
         [ 1,  1],
         [ 0,  2],
         [ 0,  0]],

        [[ 2,  2],
         [ 1,  2],
         [ 1,  1],
         [ 0,  2],
         [ 0,  0]]]])

On inspection, I see that some of the matrices are actually repeating two times, while for others it's not the case (as commented above) in output.

My questions are -

Why is this duplication happening?
How to use np.subtract.outer correctly for this example?
(or,) How to transform the matrix to the expected one? (5x2x5x2) -> (5x5x2)

Note: this is only a toy example. The original list of tuples is massive and float32.

user_901925 · Answer

No need to use outer. Just use broadcasting:

In [5]: a[:,None,:]-a[None,:,:]
Out[5]: 
array([[[ 0,  0],
        [-1,  0],
        [-1, -1],
        [-2,  0],
        [-2, -2]],

       [[ 1,  0],
        [ 0,  0],
        [ 0, -1],
        [-1,  0],
        [-1, -2]],

       [[ 1,  1],
        [ 0,  1],
        [ 0,  0],
        [-1,  1],
        [-1, -1]],

       [[ 2,  0],
        [ 1,  0],
        [ 1, -1],
        [ 0,  0],
        [ 0, -2]],

       [[ 2,  2],
        [ 1,  2],
        [ 1,  1],
        [ 0,  2],
        [ 0,  0]]])

a is (5,2), with the None extension the 2 terms are (5,1,2) and (1,5,2) which together form a (5,5,2) array.

outer can be used via:

np.subtract.outer(a,a)[:,np.arange(2),:,np.arange(2)].transpose(1,2,0)

which removes the duplicates from (5,2,5,2) and reorders the axes. But it is 3x slower.

user_652157 · Answer

Looking at the documentation for ufunc.outer (https://numpy.org/doc/stable/reference/generated/numpy.ufunc.outer.html), it appears np.subtract.outer will compute the outer product twice on the diagonal i=j. To me it doesn't look like there's a way to choose just one or the other. I recommend using the broadcasting approach, which lets you decide which dimensions you want to compute the outer product on. a[:,None]-a[None,:] gives you the same result as your itertools example.

user_10413145 · Answer

numpy.ufunc.outer works as expected. What is happening instead, that you are parsing a matrix to itertools.product and thereby lose a factor of 2 (your pairs) when calculating your loop.

numpy.ufunc.outer(A,B): flattens your input to 1D-arrays. The output would be of length len(A) x len(B). In your case that is: (5x2) x (5x2) = 100. Which is exactly the result you see.
itertools.product(A,B): will go only into the first depth of you objects and apply the product. As they consist of tuples with length two, your call to p[0]-p[1] is a subtraction of two tuples. The product itself hence is missing of 2 for the possible combinations of all elements.

Applying flatten() to your input of intertools.product will give the same result as ufunc.outer

a = np.array([(0,0), (1,1)]).flatten()
b = np.array([(2,2), (3,3)]).flatten()
c = np.subtract.outer(a, b).flatten()
d = np.asarray([p[0]-p[1] for p in product(a,b)]).flatten()

all(c==d)
> true

EDIT: To answer you questions more specifically:

The duplication occurs as numpy.ufunc.outer flattens each input.
if you insist using outer you might split your inputs, operate outer twice and merge them back together:

a = np.array([(0,0), (1,0), (1,1), (2,0), (2,2)])
b = np.asarray([p[0]-p[1] for p in product(a,a)])
c = np.subtract.outer(a[:,0], a[:,0])
d = np.subtract.outer(a[:,1], a[:,1])
e = np.vstack((c.flatten(),d.flatten())).T
np.all(b==e)
> true

Taking an outer subtraction between a list of tuples using Numpy

3 Answers