6

I have a numpy array A, which has shape (10,).

I also have, as of this moment, a numpy array B with shape (10,3,5). I want to do a multiplication between these two to get C such that C[0,:,:]=A[0]*B[0,:,:], C[1]=A[1]*B[1,:,:], etc.

I do not want to work this out with loops, one reason being the aesthetics of the thing, the other being that this code needs to be very generic. I want the user to be able to input pretty much any B of any shape as long as the leading dimension is 10. For instance, I want the user to be able to also put in a B of shape (10,4).

So: How can I implement this multiplication using numpy? Thanks.

ADDENDUM: Have been asked for example. Will go smaller. Let's say A is the numpy array [1,2,3] and B is the numpy array [[1,2],[4,5],[7,8]]. I want the multiplication of the two to result in [[1,2],[8,10],[21,24]]. ...

>>> a
array([1, 2, 3])
>>> b
array([[1, 2],
       [4, 5],
       [7, 8]])
>>> #result
>>> c
array([[ 1,  2],
       [ 8, 10],
       [21, 24]])
>>>
Improve this question
3
  • 1
    Please include, small, example arrays and the expected output. Commented Aug 24, 2015 at 18:37
  • If B is (3,5,10), the A*B works. numpy automatically adds dimensions at the start as needed (MATLAB adds them at the end). Commented Aug 24, 2015 at 18:47
  • As noted above, broadcasting pairs off the right-most dimensions of the arrays. This means that (A*B.T).T would match up the length 10s for multiplication and would be a general solution, but I think the einsum approach as suggested by DSM is arguably nicer. Commented Aug 24, 2015 at 18:53

2 Answers 2

Reset to default
4

You could use None (or np.newaxis) to expand A to match B:

>>> A = np.arange(10)
>>> B = np.random.random((10,3,5))
>>> C0 = np.array([A[i]*B[i,:,:] for i in range(len(A))])
>>> C1 = A[:,None,None] * B
>>> np.allclose(C0, C1)
True

But this will only work for the 2 case. Borrowing from @ajcr, with enough transposes we can get implicit broadcasting to work for the general case:

>>> C3 = (A * B.T).T
>>> np.allclose(C0, C3)
True

Alternatively, you could use einsum to provide the generality. In retrospect it's probably overkill here compared with the transpose route, but it's handy when the multiplications are more complicated.

>>> C2 = np.einsum('i,i...->i...', A, B)
>>> np.allclose(C0, C2)
True

and

>>> B = np.random.random((10,4))
>>> D0 = np.array([A[i]*B[i,:] for i in range(len(A))])
>>> D2 = np.einsum('i,i...->i...', A, B)
>>> np.allclose(D0, D2)
True
Improve this answer
Sign up to request clarification or add additional context in comments.

5 Comments

Whoa, some sophisticated stuff, IMHO! Actually, your first one works nicely with an edit: >>> np.array([A[i]*B[i] for i in range(len(A))]) See any problem with that? Thanks!
@bob.sacamento: no, that should work-- I was just copying your brute-force example. But since that creates the array from a slow listcomp, it will be slow for large shapes, and should only be used for testing purposes.
It's boring and inelegant, but a[(slice(None),) + (None,) * (b.ndim - 1)] * b also gets shit done...
@DSM Finally got around to trying your einsum method. I did a direct copy and paste and got the error "ValueError: operand 0 did not have enough dimensions to match the broadcasting, and couldn't be extended because einstein sum subscripts were specified at both the start and end" Any ideas? Also, know any good tutorials on einsum. I find the document page pretty opaque. Thanks.
@DSM If you're interested, I also copied and pasted an example from the einsum document -- "np.einsum('...j,j', a, b)" -- and got the same error. Any idea what this might be about? Do I have maybe an older version of scipy/numpy that doesn't handle this sort of thing?
2

Although I like the einsum notation, I'll add a little variety to the mix ....

You can add enough extra dimensions to a so that it will broadcast across b.

>>> a.shape
(3,)
>>> b.shape
(3,2)

b has more dimensions than a

extra_dims = b.ndim - a.ndim

Add the extra dimension(s) to a

new_shape = a.shape + (1,)*extra_dims    # (3,1)
new_a = a.reshape(new_shape)

Multiply

new_a * b

As a function:

def f(a, b):
    '''Product across the first dimension of b.

    Assumes a is 1-dimensional.
    Raises AssertionError if a.ndim > b.ndim or
     - the first dimensions are different
    '''
    assert a.shape[0] == b.shape[0], 'First dimension is different'
    assert b.ndim >= a.ndim, 'a has more dimensions than b'

    # add extra dimensions so that a will broadcast
    extra_dims = b.ndim - a.ndim
    newshape = a.shape + (1,)*extra_dims
    new_a = a.reshape(newshape)

    return new_a * b
Improve this answer

Comments

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.