Making the nested loop run faster, e.g. by vectorization in Python

Viewed 241

I have a 2 * N integer array ids representing intervals, where N is about a million. It looks like this

 0 2 1 ...
 3 4 3 ...

The ints in the arrays can be 0, 1, ... , M-1, where M <= 2N - 1. (Detail: if M = 2N, then the ints span all the 2N integers; if M < 2N, then there are some integers that have the same values.)

I need to calculate a kind of inverse map from ids. What I called "inverse map" is to see ids as intervals and capture the relation from their inner points with their indices.

Intuition Intuitively,

 0 2 1 
 3 4 3 

can be seen as

0 -> 0, 1, 2
1 -> 2, 3
2 -> 1, 2 

where the right-hand-side endpoints are excluded for my problem. The "inverse" map would be

0 -> 0
1 -> 0, 2
2 -> 0, 1, 2
3 -> 1

Code I have a piece of Python code that attempts to calculate the inverse map in a dictionary inv below:

  for i in range(ids.shape[1]):
    for j in range(ids[0][i], ids[1][i]):
        inv[j].append(i)

where each inv[j] is an array-like data initialized as empty before the nested loop. Currently I use python's built-in arrays to initialize it.

  for i in range(M):  inv[i]=array.array('I')

Question The nested loop above works like a mess. In my problem setting (in image processing), my first loop has a million iterations; second one about 3000 iterations. Not only it takes much memory (because inv is huge), it is also slow. I would like to focus on speed in this question. How can I accelerate this nested loop above, e.g. with vectorization?

3 Answers

avoid for loop, just a pandas sample

import numpy as np
import pandas as pd

df = pd.DataFrame({
    "A": np.random.randint(0, 100, 100000),
    "B": np.random.randint(0, 100, 100000)
})

df.groupby("B")["A"].agg(list)

You could try the below option, in which, your outer loop is hidden away within numpy's C-language implementation of apply_along_axis(). Not sure about about performance benefit, only a test at a decent scale can tell (especially as there's some initial overhead involved in converting lists to numpy arrays):

import numpy as np
import array

ids = [[0,2,1],[3,4,3]]
ids_arr = np.array(ids)   # Convert to numpy array. Expensive operation?

range_index = 0   # Initialize. To be bumped up by each invocation of my_func()

inv = {}
for i in range(np.max(ids_arr)):
    inv[i] = array.array('I')

def my_func(my_slice):
    global range_index
    for i in range(my_slice[0], my_slice[1]):
        inv[i].append(range_index)
    range_index += 1

np.apply_along_axis (my_func,0,ids_arr)
print (inv)

Output:

{0: array('I', [0]), 1: array('I', [0, 2]), 2: array('I', [0, 1, 2]), 3: array('I', [1])}

Edit:

I feel that using a dictionary might not be a good idea here. I suspect that in this particular context, dictionary-indexing might actually be slower than numpy array indexing. Use the below lines to create and initialize inv as a numpy array of Python arrays. The rest of the code can remain as-is:

inv_len = np.max(ids_arr)
inv = np.empty(shape=(inv_len,), dtype=array.array)
for i in range(inv_len):
    inv[i] = array.array('I')

(Note: This assumes that your application isn't doing dict-specific stuff on inv, such as inv.items() or inv.keys(). If that's the case, however, you might need an extra step to convert the numpy array into a dict)

Since the order of N is large, I've come up with what seems like a practical approach; let me know if there are any flaws.

  1. For the ith interval as [x,y], store it as [x,y,i]. Sort the arrays based on their start and end times. This should take O(NlogN) time.
  2. Create a frequency array freq[2*N+1]. For each interval, update the frequency using the concept of range update in O(1) per update. Generating the frequencies gets done in O(N).
  3. Determine a threshold, based on your data. According to that value, the elements can be specified as either sparse or frequent. For sparse elements, do nothing. For frequent elements only, store the intervals in which they occur.
  4. During lookup, if there is a frequent element, you can directly access the pre-computed lists. If the element is a sparse one, you can search the intervals in O(logN) time, since the intervals are sorted and there indexes were appended in step 1.

This seems like a practical approach to me, rest depends on your usage. Like the amortized time complexity you need per query and so on.

Related