7

Suppose i have:

x1 = [1, 3, 2, 4]

and:

x2 = [0, 1, 1, 0]

with the same shape

now i want to "put x2 ontop of x1" and sum up all the numbers of x1 corresponding to the numbers of x2

so the end result is:

end = [1+4 ,3+2]  # end[0] is the sum of all numbers of x1 where a 0 was in x2

this is a naive implementation using list to further clarify the question

store_0 = 0
store_1 = 0
x1 = [1, 3, 4, 2]
x2 = [0, 1, 1, 0]
for value_x1 ,value_x2 in zip(x1 ,x2):
    if value_x2 == 0:
        store_0 += value_x1
    elif value_x2 == 1:
        store_1 += value_x1

so my question: is there is a way to implement this in numpy without using loops or in general just faster?

6
  • 2
    Is it always just a few values? The numpy expression x2==1 returns a set of true/false values that can be used to filter other operations. So, x1[x2==0].sum() and x1[x2==1].sum() do the two operations you have there. Apr 26, 2021 at 18:59
  • thanks but the solution needs to be ableto handle larger arrays with more values Apr 26, 2021 at 19:00
  • Not sure why you didn't take @TimRoberts solution. I just tested with 10,000 element arrays and it took less than a second on my laptop.
    – Brad Day
    Apr 26, 2021 at 19:11
  • i meant the that the range of the x2 array could have a larger range Apr 26, 2021 at 19:13
  • was my fault that i didnt say it Apr 26, 2021 at 19:13

7 Answers 7

7

In this particular example (and, in general, for unique, duplicated, and groupby kinds of operations), pandas is faster than a pure numpy solution:

A pandas way, using Series (credit: very similar to @mcsoini's answer):

def pd_group_sum(x1, x2):
    return pd.Series(x1, index=x2).groupby(x2).sum()

A pure numpy way, using np.unique and some fancy indexing:

def np_group_sum(a, groups):
    _, ix, rix = np.unique(groups, return_index=True, return_inverse=True)
    return np.where(np.arange(len(ix))[:, None] == rix, a, 0).sum(axis=1)

Note: a better pure numpy way is inspired by @Woodford's answer:

def selsum(a, g, e):
    return a[g==e].sum()

vselsum = np.vectorize(selsum, signature='(n),(n),()->()')

def np_group_sum2(a, groups):
    return vselsum(a, groups, np.unique(groups))

Yet another pure numpy way is inspired by a comment from @mapf about using argsort(). That in itself already takes 45ms, but we may try something based on np.argpartition(x2, len(x2)-1) instead, since that takes only 7.5ms by itself on the benchmark below:

def np_group_sum3(a, groups):
    ix = np.argpartition(groups, len(groups)-1)
    ends = np.nonzero(np.diff(np.r_[groups[ix], groups.max() + 1]))[0]
    return np.diff(np.r_[0, a[ix].cumsum()[ends]])

(Slightly modified) example

x1 = np.array([1, 3, 2, 4, 8])  # I added a group for sake of generality
x2 = np.array([0, 1, 1, 0, 7])

>>> pd_group_sum(x1, x2)
0    5
1    5
7    8

>>> np_group_sum(x1, x2)  # and all the np_group_sum() variants
array([5, 5, 8])

Speed

n = 1_000_000
x1 = np.random.randint(0, 20, n)
x2 = np.random.randint(0, 20, n)

%timeit pd_group_sum(x1, x2)
# 13.9 ms ± 65.6 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

%timeit np_group_sum(x1, x2)
# 171 ms ± 129 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

%timeit np_group_sum2(x1, x2)
# 66.7 ms ± 19.4 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

%timeit np_group_sum3(x1, x2)
# 25.6 ms ± 41.3 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

Going via pandas is faster, in part because of numpy issue 11136.

11
  • was very confused when pandas was 10x faster than numpy Apr 26, 2021 at 20:02
  • but pandas seems to be the best way Apr 26, 2021 at 20:02
  • yes; because of the np.unique() slowness issue, the other day I went through the pandas code for .duplicated() and family. It is quite optimized indeed. In fact, see my comment at the end of that numpy issue referenced above.
    – Pierre D
    Apr 26, 2021 at 20:08
  • @PierreD I just wanted to point out that I feel like there is a smart way hidden somwhere using sorting (e.g., x2inds = np.array(x2).argsort(); x1_sorted = np.array(x1)[x2inds] or similar) but I didn't manage to find it. Maybe you know a way. Also it would be nice if you include the other approaches mentioned here in your speed comparison! :)
    – mapf
    Apr 26, 2021 at 20:23
  • 2
    @mapf: on the "benchmark" above, just %timeit np.argsort(x2) already takes 43ms, and np.unique(x2) takes 31ms. So any approach built on those will lose against pandas. But following that line of thought, I tried a version using np.argpartition(), as that takes only 7.5ms. I updated my answer to include a way using that, but it still takes 26ms.
    – Pierre D
    Apr 26, 2021 at 21:06
5
>>> x1 = np.array([1, 3, 2, 7])
>>> x2 = np.array([0, 1, 1, 0])
>>> for index in np.unique(x2):
>>>     print(f'{index}: {x1[x2==index].sum()}')
0: 8
1: 5
>>> # or in one line
>>> [(index, x1[x2==index].sum()) for index in np.unique(x2)]
[(0, 8), (1, 5)]
5
  • Small idea for improvement: I'm not sure, but I guess set(x2) would be a bit faster than np.unique(x2).
    – mapf
    Apr 26, 2021 at 19:09
  • i forgot to say that the solution should be able to handle "store" values ,the best thing would be a function that returns an array with the corresponding values.So the x2 could also have a range from 0 to 1000 Apr 26, 2021 at 19:10
  • @mapf Could be, but depends how numpy implements unique. I generally trust numpy to be efficient.
    – Woodford
    Apr 26, 2021 at 19:11
  • 1
    @user15770670 This code will handle any size arrays and and any number of indices in x2. Not sure what else you're asking for.
    – Woodford
    Apr 26, 2021 at 19:14
  • with np.vectorize(), this is better than my initial pure-numpy answer! Kudos.
    – Pierre D
    Apr 26, 2021 at 20:29
3

Would a pandas one-liner be ok?

store_0, store_1 = pd.DataFrame({"x1": x1, "x2": x2}).groupby("x2").x1.sum()

Or as a dictionary, for arbitrarily many values in x2:

pd.DataFrame({"x1": x1, "x2": x2}).groupby("x2").x1.sum().to_dict()

Output:

{0: 5, 1: 5}
3
  • i forgot to say that the solution should be able to handle "store" values ,the best thing would be a function that returns an array with the corresponding values.So the x2 could also have a range from 0 to 1000 Apr 26, 2021 at 19:11
  • didnt expect that but pandas outperformed numpy by a faktor of 8x!!! Apr 26, 2021 at 19:45
  • one of the fastest ways and faster than pure numpy. Good stuff.
    – Pierre D
    Apr 26, 2021 at 20:30
2

using compress

from itertools import compress
result = [sum(compress(x1,x2)),sum(compress(x1, (map(lambda x: not x,x2))))]
1

This extends your loop into a larger number of values. I can't think of a numpy one-liner to do this.

sums = [0] * 10000
for vx1,vx2 in zip(x1,x2):
    sums[vx2] += vx1
2
  • i forgot to say that the solution should be able to handle "store" values ,the best thing would be a function that returns an array with the corresponding values.So the x2 could also have a range from 0 to 1000 Apr 26, 2021 at 19:10
  • im dumb :( i think this will do it! Apr 26, 2021 at 19:17
1

By casting the second list as a Boolean array, you can use it to index the first one:

import numpy as np

x1 = np.array([1, 3, 2, 4])
x2 = np.array([0, 1, 1, 0], dtype=bool)

end = [np.sum(x1[~x2]), np.sum(x1[x2])]
end
[5, 5]

Edit: If x2 can have values larger than 1, you could use a list comprehension:

x1 = np.array([1, 3, 2, 4])
x2 = np.array([0, 1, 1, 0])

end = [np.sum(x1[x2 == i]) for i in range(max(x2) + 1)]
1

This extends the solution Tim Roberts suggested at the begining but will account for X2 having multiple values i.e Non binary. Here those values are strictly adjacent because the for loop uses the range of rng but it could be extended so that x2 has values that are not adjacent e.g [0 2 2 2 1 4] <- no 3's whereas randint used for this example will return a vector something like [0 1 1 3 4 2].

import numpy as np
rng = 5 # Range of values for x2 i.e [0 1 2 3 4]
x1 = np.random.randint(20, size=10000) #random vector of size 10k
x2 = np.random.randint(5, size=10000) # inexing vector size 10k with range (0-4)


store = []
for i in range(rng): # loop and append to list
    store.append(x1[x2==i].sum()) 
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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