Building 3D arrays in Python to replace loops for optimization

Question

I'm trying to better understand python optimization so this is a dummy case, but hopefully outlines my idea...

Say I have a function which takes two variables:

def func(param1, param2):
    return some_func(param1) + some_const*(param2/2)

and I have arrays for param1 and param2 (of different lengths), at which I want the function to be evaluated, (some_func is an arbitrary function of param1) e.g.

param1 = np.array((1,2,3,4,5))
param2 = np.array((5,2,3,1,9, 9, 10))

I can evaluate over all parameter space by doing:

result = []
for p in param1:
    result.append(func(p, param2))
result = np.asarray(result)

However, loops in Python are slower than array operations. Therefore, I wonder is there a way to achieve a 3D array which contains the results of func for all values in both param1 and param2 arrays?

hpaulj · Accepted Answer · 2014-08-03 05:40:26Z

3

Your example:

In [198]: result=[]
In [199]: for p in param1:
   .....:     result.append(p+3*(param2/2))
In [200]: result=np.array(result)

Same result using broadcasting (and np.newaxis):

In [197]: param1[:,None] + 3*(param2[None,:]/2)
Out[197]: 
array([[ 7,  4,  4,  1, 13, 13, 16],
       [ 8,  5,  5,  2, 14, 14, 17],
       [ 9,  6,  6,  3, 15, 15, 18],
       [10,  7,  7,  4, 16, 16, 19],
       [11,  8,  8,  5, 17, 17, 20]])

Details of some_func would determine whether you use some_func(param1[:,None]) or some_func(param1)[:,None].

answered Aug 3, 2014 at 5:40

hpaulj

233k14 gold badges260 silver badges392 bronze badges

Sign up to request clarification or add additional context in comments.

Comments

Roland Smith · Accepted Answer · 2014-08-03 10:11:32Z

3

Original answer for some_func(param1) x param2

Write the some_func in such a way that it can accept and return numpy arrays. Then use;

numpy.outer(some_func(param1), param2)

This works because in your example, both param1 and param2 are vectors (1D arrays), so you can use outer and the result will be a 2D array, not 3D.

Edit

As long as the operation you want to do is a universal function ("ufunc"), you can use its outer method;

In [1]: import numpy as np

In [2]: a = np.arange(10)

In [3]: a
Out[3]: array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

In [4]: b = np.arange(15)

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

edited Aug 3, 2014 at 10:11

answered Aug 2, 2014 at 22:59

Roland Smith

43.8k3 gold badges69 silver badges98 bronze badges

5 Comments

IanRoberts Over a year ago

Thanks for your response. This would work for the specific case of multiplication as I originally put in my question. In fact, I'm interested in this more generally when the function isn't a simple multiplication (e.g. return some_func(param1) + some_const*(param2/2)) as in this case I imagine np.outer would no longer work. Thanks

Eelco Hoogendoorn Over a year ago

indeed. outer does multiplication, whereas addition is called for

Roland Smith Over a year ago

@IanRoberts All universal functions support outer. See improved answer. All standard arithmatic operations are available as ufuncs.

IanRoberts Over a year ago

Thanks for clarifying, is this essentially the same as the procedure outline in hpaulji's answer? And is there any difference in terms of memory usage of speed by using the numpy outer function rather than broadcasting and np.newaxis?

Roland Smith Over a year ago

@IanRoberts I'm not sure. It could be that usinf ufunc,outer directly is faster because it looks like it uses less operations. But you'd have to test it.

Eelco Hoogendoorn · Accepted Answer · 2014-08-03 07:11:22Z

2

You can use numpy broadcasting, as hpaulj suggested. This is generally the first line of attack for vectorizing numpy operations. The essence of using broadcasting, is that in doing so, you move your loops from the slow python level down to the C level.

Alternatively, you may choose to take a look at numba, and numexpr, which accomplish vectorization in a more flexible manner

answered Aug 3, 2014 at 7:11

Eelco Hoogendoorn

10.8k1 gold badge46 silver badges43 bronze badges

Collectives™ on Stack Overflow

Building 3D arrays in Python to replace loops for optimization

3 Answers 3

Comments

5 Comments

Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

3 Answers 3

Comments

5 Comments

Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related