What makes the different performances between numpy.sum and numpy.cumsum?

Viewed 745

If my understanding is correct, both np.sum and np.cumsum will take O(n) time. However, when I do both operations to the same matrix sequentially (on different axis), it seems the sequence of np.sum and np.cumsum makes the overall performance quite different, though the result is the same as expected.

If I do the np.cumsum along column direction (axis=1) for each row first, and then do np.sum for all the rows (axis=0), it will take longer time.

If I do the np.sum along rows (axis=0) first, then do np.cumsum for the one dimension array, it will take shorter time.

My hypothesis is that the np.cumsum will take more time on the data allocation/manipulation since it will yield more data than np.sum, so it will take longer time if there are more operations of np.cumsum.

Here is my testing code and result

import numpy as np
import time

b = np.zeros((1000, 1000))
for i in range(1000):
    b[i] = np.array(range(1000))

time_start = time.time()
for i in range(1000):
    c = np.cumsum(b, axis=1)
    d = np.sum(c, axis=0)
time_end = time.time()
print(f"np.sum(np.cumsum(...)) time: {time_end - time_start}")


time_start = time.time()
for i in range(1000):
    c = np.cumsum(np.sum(b, axis=0))
time_end = time.time()
print(f"np.cumsum(np.sum(...)) time: {time_end - time_start}")

The result is

np.sum(np.cumsum(...)) time: 3.6612446308135986
np.cumsum(np.sum(...)) time: 0.38796162605285645
2 Answers

Look at your alternatives applied to a small array.

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

sum(cumsum(b)):

In [47]: np.cumsum(b, axis=1)         # same shape as b
Out[47]: 
array([[  0,   1,   3,   6,  10,  15],
       [  6,  13,  21,  30,  40,  51],
       [ 12,  25,  39,  54,  70,  87],
       [ 18,  37,  57,  78, 100, 123]])
In [48]: np.sum(_, axis=0)
Out[48]: array([ 36,  76, 120, 168, 220, 276])

sum first:

In [49]: np.sum(b, axis=0)            # much reduced array
Out[49]: array([36, 40, 44, 48, 52, 56])
In [50]: np.cumsum(_)
Out[50]: array([ 36,  76, 120, 168, 220, 276])

While the first step operates on the same number of elements, the whole b, the second step works with a much reduced array when sum is done first.

Actually for this small sample the time differences are minor:

In [57]: timeit np.sum(np.cumsum(b,axis=1), axis=0)
17.8 µs ± 29.5 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
In [58]: timeit np.cumsum(np.sum(b, axis=0))
16.1 µs ± 30.5 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

But for a bigger array, the time savings due to the reduction becomes significant:

In [59]: b = np.random.randint(0,1000,(1000,1000))
In [60]: timeit np.sum(np.cumsum(b,axis=1), axis=0)
5.11 ms ± 9.99 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
In [61]: timeit np.cumsum(np.sum(b, axis=0))
1.64 ms ± 649 ns per loop (mean ± std. dev. of 7 runs, 1000 loops each)

Looking at the first step by itself:

In [62]: timeit np.cumsum(b,axis=1)
3.73 ms ± 12.6 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
In [63]: timeit np.sum(b, axis=0)
1.62 ms ± 1.7 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

cumsum, because it returns a full b.shape array is slower.

O(n) tells us something about how an operation scales, such as when operating on small arrays versus big ones. But it can't be used to compare one operation to another. In addition with numpy, the core calculations are done rapidly in compiled code, but memory management and other setup issues may swamp that core time.

edit

np.sum is essentially np.add.reduce:

In [79]: timeit np.add.reduce(b,axis=0)
1.59 ms ± 4.73 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
In [80]: timeit np.sum(b,axis=0)
1.64 ms ± 5.48 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

And np.cumsum is np.add.accumulate:

In [81]: timeit np.cumsum(b,axis=0)
10.2 ms ± 357 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
In [82]: timeit np.add.accumulate(b,axis=0)
10 ms ± 14.6 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

If my understanding is correct, both np.sum and np.cumsum will take O(n) time.

Correct (assuming the input is the same)!

it seems the sequence of np.sum and np.cumsum makes the overall performance quite different

np.cumsum is more expensive than np.sum on the same input since np.cumsum needs to allocate and write an output array (which can be pretty big). Moreover, assuming floating-point operations are associative, the implementation of np.sum can be easily optimized while it is quite hard to optimize np.cumsum.

My hypothesis is that the np.cumsum will take more time on the data allocation/manipulation since it will yield more data than np.sum, so it will take longer time if there are more operations of np.cumsum.

The thing is that the first implementation have much more work to do than the second. This is mainly why it is slower, especially because of the memory reads/writes. Indeed, np.cumsum produces a big 2D array in the first implementation that must be computed by np.sum. Writing and reading this temporary array is expensive. In the example, the first implementation require to move 14.9 GiB from/to the memory hierarchy (likely in RAM) just for the temporary array. Note that building a temporary array also make the computation less cache-friendly as it requires more space and thus b and c may not both fit in the cache any more (b and c takes 8 MiB each on my machine and my CPU have a L3 cache of 9 MiB which mean that not both can fit in the cache) as opposed to the second implementation. The throughput of the RAM is often much smaller than the one of the CPU caches.

Note that the second implementation also produces temporary arrays. This version should allocate as much temporary arrays as the first one. However, the second implementation produces temporary arrays that are 1000 time smaller! Thus, the call to np.cumsum is much faster in this version. On my machine, b is mostly stored in the fast L3 cache and temporary arrays in the very fast L1 cache.

Related