Vectorizing Nested Loop With Conditionals And Functions

I have the following function: def F(x): #F receives a numpy vector (x) with size (xsize*ysize) ff = np.zeros(xsize*ysize) count=0 for i in range(xsize):

Solution 1:

The first thing to note is that your function F(x) can be described as x(idx) * weight(idx) for each index, where weight is only dependent on the dimensions of x. So lets structure our code in terms of a function get_weights_for_shape so that F is fairly simple. For sake of simplicity weights will be a (xsize by size) matrix but we can let F work for flat inputs too:

defF(x, xsize=None, ysize=None):
    iflen(x.shape) == 2:
        # based on how you have put together your question this seems like the most reasonable representation.
        weights = get_weights_for_shape(*x.shape)
        return x * weights
    eliflen(x.shape) == 1and xsize * ysize == x.shape[0]:
        # single dimensional input with explicit size, use flattened weights.
        weights = get_weights_for_shape(xsize, ysize)
        return x * weights.flatten()
    else:
        raise TypeError("must take 2D input or 1d input with valid xsize and ysize")


# note that get_one_weight=function can be replaced with your actual function.defget_weights_for_shape(xsize, ysize, get_one_weight=function):
    """returns weights matrix for F for given input shape"""# will use (xsize, ysize) shape for these calculations.
    weights = np.zeros((xsize,ysize))
    #TODO: will fill in calculations herereturn weights

So first we want to run your function (which I have aliased get_one_weight inside the parameters) for each element, you have said that this function can't be vectorized so we can just use list comprehension. We want a matrix a that has the same shape (xsize,ysize) so the comprehension is a bit backwards for a nested list:

# notice that the nested list makes the loops in opposite order:# [ROW for i in Xs]#  ROW = [f() for j in Ys]a = np.array([[get_one_weight(i,j,xsize,ysize)
                    for j in range(ysize)
              ] for i in range(xsize)])

With this matrix a > xsize will give a boolean array for conditional assignment:

case1 = a > xsize
weights[case1] = a[case1]

For the other case, we use the index i and j. To vectorize a 2D index we can use np.meshgrid

[i,j] = np.meshgrid(range(xsize), range(ysize), indexing='ij')
case2 = ~case1 # could have other cases, in this case it's just the rest.
weights[case2] = i[case2] * j[case2]

return weights #that covers all the calculations

Putting it all together gives this as the fully vectorized function:

# note that get_one_weight=function can be replaced with your actual function.defget_weights_for_shape(xsize, ysize, get_one_weight=function):
    """returns weights matrix for F for given input shape"""# will use (xsize, ysize) shape for these calculations.
    weights = np.zeros((xsize,ysize))

    # notice that the nested list makes the loop order confusing:# [ROW for i in Xs]#  ROW = [f() for j in Ys]
    a = np.array([[get_one_weight(i,j,xsize,ysize)
                        for j inrange(ysize)
                  ] for i inrange(xsize)])

    case1 = (a > xsize)
    weights[case1] = a[case1]

    # meshgrid lets us use indices i and j as vectorized matrices.
    [i,j] = np.meshgrid(range(xsize), range(ysize), indexing='ij')
    case2 = ~case1
    weights[case2] = i[case2] * j[case2]
    #could have more than 2 cases if applicable.return weights

And that covers most of it. For your specific case since this heavy calculation only relies on the shape of the input, if you expect to call this function repeatedly with similarly sized input you can cache all previously calculated weights:

defget_weights_for_shape(xsize, ysize, _cached_weights={}):
    if (xsize, ysize) notin _cached_weights:
        #assume we added an underscore to real function written above
        _cached_weights[xsize,ysize] = _get_weights_for_shape(xsize, ysize)
    return _cached_weights[xsize,ysize]

As far as I can tell this seems like the most optimized you are going to get. The only improvements would be to vectorize function (even if that means just calling it in several threads in parallel) or possibly if .flatten() makes an expensive copy that could be improved but I'm not totally sure how.