Take the 2-minute tour ×
Code Review Stack Exchange is a question and answer site for peer programmer code reviews. It's 100% free, no registration required.

I need to speed up the function below:

import numpy as np
import itertools

def combcol(myarr):
    ndims = myarr.shape[0]
    solutions = []
    for idx1, idx2, idx3, idx4, idx5, idx6 in itertools.combinations(np.arange(ndims), 6):
        c1, c2, c3, c4, c5, c6 = myarr[idx1,1], myarr[idx2,2], myarr[idx3,1], myarr[idx4,2], myarr[idx5,1], myarr[idx6,2]
        if c1-c2>0 and c2-c3<0 and c3-c4>0 and c4-c5<0 and c5-c6>0 :
            solutions.append(((idx1, idx2, idx3, idx4, idx5, idx6),(c1, c2, c3, c4, c5, c6)))
    return solutions

X = np.random.random((20, 10))  
Y = np.random.random((40, 10))  

if __name__=='__main__':
    from timeit import Timer
    t = Timer(lambda : combcol(X))
    t1 = Timer(lambda : combcol(Y))
    print('t : ',t.timeit(number=1),'t1 : ',t1.timeit(number=1))

Results:

t :  0.6165180211451455 t1 :  64.49216925614847 

The algorithm is too slow for my standard use (myarr.shape[0] = 500). Is there a NumPy way to decrease execution time of this function (without wasting too much memory)? Is it possible to implement the problem in Cython?

share|improve this question
2  
What are you trying to achieve here? –  Gareth Rees Dec 29 '13 at 20:08

1 Answer 1

up vote 7 down vote accepted

Your algorithm is of time complexity O(n^6) in ndims. If your current code needs 65 s for ndims = 40, then one could estimate the time for ndims = 500 to be roughly 8 years. Even an improvement of the most inner block's execution time by a factor of 1000 would yield an execution time of about 3 days. I guess this is not good enough.

So you should better try to improve the algorithm itself. However I do not really understand the problem the function is trying to solve and I couldn't come up with a better solution on the spot.

Also I haven't looked at numpy approaches, because you seem to have no vectorizable calculations in there.

Nonetheless I tried to make improvements to the code as is. I came up with the following:

  1. You use only the a small fraction of myarr, so extract these parts first

    myarr1 = myarr[:,1]
    myarr2 = myarr[:,2]
    
  2. Since you do many index lookups on these arrays turn them into lists. Numpy is not good at single index lookups.

    myarr1 = myarr[:,1].tolist()
    myarr2 = myarr[:,2].tolist()
    
  3. Instead of substracting and comparing to 0, compare the variables directly saving one operation:

    if c1>c2 and c2<c3 and c3>c4 and c4<c5 and c5>c6:
    
  4. Instead of using and, you can also connect the comparisons all together, potentially saving the double variable lookups:

     if c1 > c2 < c3 > c4 < c5 > c6:
    
  5. Instead of assigning the array elements to c1, c2, etc., use them directly. Since later comparisons are not executed if earlier one already return False this might save array lookups:

    if myarr1[idx1] > myarr2[idx2] < myarr1[idx3] > myarr2[idx4] < myarr1[idx5] > myarr2[idx6]:
    

All in all I get the following code:

import numpy as np
import itertools

def combcol(myarr):
    ndims = myarr.shape[0]
    solutions = []
    for idx1, idx2, idx3, idx4, idx5, idx6 in itertools.combinations(np.arange(ndims), 6):
        c1, c2, c3, c4, c5, c6 = myarr[idx1,1], myarr[idx2,2], myarr[idx3,1], myarr[idx4,2], myarr[idx5,1], myarr[idx6,2]
        if c1-c2>0 and c2-c3<0 and c3-c4>0 and c4-c5<0 and c5-c6>0 :
            solutions.append(((idx1, idx2, idx3, idx4, idx5, idx6),(c1, c2, c3, c4, c5, c6)))
    return solutions

def combcol2(myarr):
    ndims = myarr.shape[0]
    myarr1 = myarr[:,1].tolist()
    myarr2 = myarr[:,2].tolist()
    solutions = []
    for idx1, idx2, idx3, idx4, idx5, idx6 in itertools.combinations(range(ndims), 6):
        if myarr1[idx1] > myarr2[idx2] < myarr1[idx3] > myarr2[idx4] < myarr1[idx5] > myarr2[idx6]:
            solutions.append(((idx1, idx2, idx3, idx4, idx5, idx6),(myarr1[idx1], myarr2[idx2], myarr1[idx3], myarr2[idx4], myarr1[idx5], myarr2[idx6])))
    return solutions

X = np.random.random((40, 10))
Y = [np.random.random((10, 10)) for i in range(100)]

if __name__=='__main__':
    from timeit import Timer
    for test_case in Y: assert combcol(test_case) == combcol2(test_case)
    t = Timer(lambda : combcol(X))
    t1 = Timer(lambda : combcol2(X))
    print('t : ',t.timeit(number=1),'t1 : ',t1.timeit(number=1))

With results like these:

t :  12.221544027328491 t1 :  0.8104081153869629
t :  10.69440507888794 t1 :  0.6592199802398682
t :  12.91178011894226 t1 :  0.8727250099182129
t :  11.21804690361023 t1 :  0.6851639747619629

The improvement is not even close to what you need. Also I tried to run the algorithm with ndims = 80 and got memory problems. The reason is the number of solutions: There were 18679356 solutions taking up 5 GB of RAM. With other values for ndims I found that around 5% - 25% of all possible solutions are valid! So it seems you need to change the way you output your solutions if you want to get beyond ndims = 100 on any realistic machine.

Edit: I just notice that 5 GB is a bit too much for the number of solutions. However this was the memory allocated by the python process. A run with ndims = 100 was unsuccessful due to lack of available memory.

Edit2: There is actually something far more efficient you can do, by interweaving your solution checks with the combination generation (abandoning itertools.combinations). This doesn't look so nice however:

import numpy as np
import itertools


def combcol(myarr):
    ndims = myarr.shape[0]
    solutions = []
    for idx1, idx2, idx3, idx4, idx5, idx6 in itertools.combinations(np.arange(ndims), 6):
        c1, c2, c3, c4, c5, c6 = myarr[idx1,1], myarr[idx2,2], myarr[idx3,1], myarr[idx4,2], myarr[idx5,1], myarr[idx6,2]
        if c1-c2>0 and c2-c3<0 and c3-c4>0 and c4-c5<0 and c5-c6>0 :
            solutions.append(((idx1, idx2, idx3, idx4, idx5, idx6),(c1, c2, c3, c4, c5, c6)))
    return solutions


def combcol2(myarr):
    ndims = myarr.shape[0]
    myarr1 = myarr[:,1].tolist()
    myarr2 = myarr[:,2].tolist()
    solutions = []

    for idx1 in range(ndims):
        c1 = myarr1[idx1]
        for idx2 in range(idx1+1, ndims):
            c2 = myarr2[idx2]
            if c1 > c2:
                for idx3 in range(idx2+1, ndims):
                    c3 = myarr1[idx3]
                    if c2 < c3:
                        for idx4 in range(idx3+1, ndims):
                            c4 = myarr2[idx4]
                            if c3 > c4:
                                for idx5 in range(idx4+1, ndims):
                                    c5 = myarr1[idx5]
                                    if c4 < c5:
                                        for idx6 in range(idx5+1, ndims):
                                            c6 = myarr2[idx6]
                                            if c5 > c6:
                                                solutions.append(((idx1, idx2, idx3, idx4, idx5, idx6),(c1, c2, c3, c4, c5, c6)))
    return solutions


X = np.random.random((40, 10))
Y = [np.random.random((10, 10)) for i in range(100)]

if __name__=='__main__':
    from timeit import Timer
    for test_case in Y: assert combcol(test_case) == combcol2(test_case)
    t = Timer(lambda : combcol(X))
    t1 = Timer(lambda : combcol2(X))
    print('t : ',t.timeit(number=1),'t1 : ',t1.timeit(number=1))

The result:

t :  11.706523180007935 t1 :  0.15908312797546387
t :  12.531341075897217 t1 :  0.16385602951049805
t :  11.581043004989624 t1 :  0.14388418197631836

So with that we already have nearly a factor of 100x speed-up.

share|improve this answer
    
Thank you very much. –  querzy Dec 30 '13 at 13:12

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.