Avoiding inner for loops??

[Martin De Kauwe]

Hi, 

I need to avoid (at least) two inner for loops in what I am trying to do otherwise my processing takes forever. What is the best way to transfer what I am doing into a more "numpy way"? Essentially I am trying to call a model again for various different parameter combinations. The example is fictional, by the grid_size would ideally grow > 500 and by doing so the processing speed becomes very slow the way I have set things up.. 

thanks. 

example. 


import numpy as np 

def fake_model(data1, data2, p1, p2, p3): 
    """ complete nonsense """ 
    return data1 + data2 * p1 * p2 * p3 

data1 = np.random.rand(10) # the size of this arrays varies might be 10 might be 15 etc 
data2 = np.random.rand(10) # the size of this arrays varies might be 10 might be 15 etc 
obs = np.random.rand(10) # the size of this arrays varies might be 10 might be 15 etc 

grid_size = 10 # Ideally this would be a large number 
param1 = np.linspace(5.0, 350, grid_size) 
param2 = np.linspace(5.0, 550, grid_size) 
param3 = np.linspace(1E-8, 10.5, grid_size) 
ss = np.zeros(0) 

for p1 in param1: 
    for p2 in param2: 
        for p3 in param3: 
            ans = fake_model(data1, data2, p1, p2, p3) 
            
            ss = np.append(ss, np.sum(obs - ans)**2) 
            print np.sum(obs - ans)**2

[Martin De Kauwe]

Perhaps simplifying, in 2D this is what I want if using loops

def fake_model(data1, data2, p1, p2): 

    """ complete nonsense """ 

    return data1 + data2 * p1 * p2

grid_size = 5

nobs = 5

obs = np.zeros(nobs)

data1 = np.arange(nobs)

data2 = np.arange(nobs)

a = np.arange(grid_size)

b = np.arange(grid_size)

c = np.arange(grid_size)

ss = np.zeros(0) 

for p1 in a: 

    for p2 in b: 

        ans = fake_model(data1, data2, p1, p2) 

        #ss = np.append(ss, np.sum(obs - ans)**2) 

        print ans

which would produce

[0 1 2 3 4]

[0 1 2 3 4]

[0 1 2 3 4]

[0 1 2 3 4]

[0 1 2 3 4]

[0 1 2 3 4]

[0 2 4 6 8]

.

snip

.

[0 1 2 3 4]

[ 0  5 10 15 20]

[ 0  9 18 27 36]

[ 0 13 26 39 52]

[ 0 17 34 51 68]

And so I figured something like...

a = np.ones((grid_size,grid_size)) * np.arange(grid_size)[None,:] 

b = np.ones((grid_size,grid_size)) * np.arange(grid_size)[:,None] 

ans = fake_model(data1, data2, a, b)

Although this doesn't seem to work, but I think this might be along the right lines? This produces

[[  0.   1.   2.   3.   4.]

 [  0.   2.   6.  12.  20.]

 [  0.   3.  10.  21.  36.]

 [  0.   4.  14.  30.  52.]

 [  0.   5.  18.  39.  68.]]

[Martin De Kauwe]

Or perhaps an easier solution would be to build the sampling grid first and then use a single for loop to run the model?

import itertools

import numpy as np

grid_size = 500

a = np.linspace(5.0, 350, grid_size) 

b = np.linspace(5.0, 550, grid_size) 

c = np.linspace(1E-8, 10.5, grid_size)

x = []

for (i,j,k) in itertools.product(a, b, c):

    x.extend((i,j,k))

This would achieve what I want but is again very slow, so is there a way to jump over the need for the two inner loops?

thanks

[Luuk van der Velden]

Consider using broadcasting to three 3D parameters matrices (A,B,C), then creating a 'ufunc' that takes three parameter at a time (a,b,c) for every identical position in the three arrays (A,B,C). So giving the broadcasting arrays as input to a ufunction which maps the model function on the three paramater arrays.

greets,

Luuk