How to apply cummulative sum in Numpy in slices with condition to previous value?

Question

I have a vector with signals with values of 1 or -1. I want to have a second vector that computes the cumulative sum of consecutive signals with the same value and restarts the cumulative sum every time the signal change. Here is an example:

signal  = [1  1  1 -1 -1 -1 -1]

cum_sum = [1  2  3 -1 -2 -3 -4]

I have large data to be computed and want to do it as efficiently as possible. My code right now does the job but it takes time and is not taking advantage of numpy efficiency:

import numpy as np

# Signal values to be analyzed
signal = np.array([1,1,1,-1,-1,-1,-1], dtype=int)

# Vector with previous value of signal
signal_prev = signal[:-1]
signal_prev = np.pad(signal_prev,(1,0), mode='constant', constant_values=(0))

#Array with signal values in first column and previous values in second column 
arr = np.array([signal,signal_prev], dtype=int)
arr = np.transpose(arr)

print(arr)
""" Array with signal values and previous values
[[ 1  0]
 [ 1  1]
 [ 1  1]
 [-1  1]
 [-1 -1]
 [-1 -1]
 [-1 -1]]
"""

#create an empty array to append cumulative sum
signal_sum = np.array([], dtype=int)

# compute the cumulative sum iterating row by row
for x in arr:
    if np.sign(x[0]*x[1]) > 0:
        signal_sum = np.append(signal_sum, signal_sum[-1] + x[1])
    else:
        signal_sum= np.append(signal_sum, x[0])

arr_sum = np.array([signal, signal_sum])
arr_sum = np.transpose(arr_sum)
print(arr_sum)
""" Array with signal values and cumulative sum restarted with signal change
[[ 1  1]
 [ 1  2]
 [ 1  3]
 [-1 -1]
 [-1 -2]
 [-1 -3]
 [-1 -4]]
"""

I believe that this calculation can be done more efficiently using numpy functions or using lambda functions. I'm not a programmer, and I'm new to Python. I would like to know if this could be done faster.

user_758174 · Answer

For a fast, fully vectorized way (no loops), you can use a regular np.cumsum(), but on a copy of your array where you subtract the previous group sum at the start of each group:

def group_cumsum(s):
    # make a copy and ensure np.array (in case list was given)
    s = np.array(s).copy()
    idx = np.nonzero(np.diff(s))[0]  # last of each group
    off = np.diff(np.concatenate(([0], np.cumsum(s)[idx])))
    s[idx + 1] -= off
    return np.cumsum(s)

Example:

print(group_cumsum([1, 1, 1, -1, -1, -1, -1]))
# [ 1  2  3 -1 -2 -3 -4]

print(group_cumsum([1]*3 + [-1]*2 + [1]*4 + [-1]*5))
# [ 1  2  3 -1 -2  1  2  3  4 -1 -2 -3 -4 -5]

The time saving is substantial for large arrays:

no loops in the Python code, all ops are vectorized, and
it is O(n + k) for k groups in an array of size n (unlike other solutions that are O(n * k)).

Try this:

s = np.random.choice([1, -1], size=(int(1e6)))

%%timeit
group_cumsum(s)

19.1 ms ± 137 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

user_10733051 · Answer

A little hacky but I am not sure if this meets your speed requirements-

import numpy as np

signal = np.array([-1, -1, 1, -1, 1, 1, -1, -1, -1, 1, 1, 1, -1, -1, 1, -1, 1, -1, 1, -1])

normal_cumsum = np.cumsum(signal)

reset_sums = np.roll(np.where(np.diff(np.sign(signal), append=np.sign(signal[-1])) != 0, normal_cumsum, np.zeros_like(signal)), 1)
# for forward fill operation
ffill_idxs = np.hstack((np.squeeze(np.argwhere(np.diff(signal, prepend=signal[0]) != 0)), np.array([len(signal)])))

for start, end in zip(ffill_idxs[:-1], ffill_idxs[1:]):
    reset_sums[start:end] = reset_sums[start]

result = normal_cumsum - reset_sums
print(result)

[-1 -2  1 -1  1  2 -1 -2 -3  1  2  3 -1 -2  1 -1  1 -1  1 -1]

What size of lists are you working with? Maybe a comparison can be made to see the efficiency

user_7153928 · Answer

2 possible approaches. sawtooth1 uses the same idea as @sai.

import numpy as np

arr = np.array([ 1,1,1,1,1,-1,-1,-1,-1,1,1,1,1,-1,-1,-1])
arr
# Out[3]: array([ 1,  1,  1,  1,  1, -1, -1, -1, -1,  1,  1,  1,  1, -1, -1, -1])


def sawtooth0( arr ):
    cumsum_ix = np.zeros_like(arr)
    cumsum_ix[ 1: ] = (( arr[ :-1 ] - arr[ 1: ] ) != 0 ).cumsum()
    # cumsum_ix incrememnts for each sign change. 

    result = np.zeros_like( arr )
    for i in range( 0, cumsum_ix[ -1 ] + 1 ):
        # For each cumsum_ix select those items and generate the cumsums.
        result[ cumsum_ix==i ] = arr[ cumsum_ix==i ].cumsum()
    return result

sawtooth0( arr )
# Out[21]: array([ 1,  2,  3,  4,  5, -1, -2, -3, -4,  1,  2,  3,  4, -1, -2, -3])

def sawtooth1( arr ):
    cumsum_ix = np.zeros_like(arr)
    cumsum_ix[ 1: ] = (( arr[ :-1 ] - arr[ 1: ] ) != 0 ).cumsum()
    # cumsum_ix incrememnts for each sign change. 

    totals = np.zeros( cumsum_ix[ -1 ] + 1, dtype = np.int )
    for i in range( 0, cumsum_ix[ -1 ] ):
        # For each cumsum_ix select those items and generate the sum.
        totals[ i+1 ] = arr[ cumsum_ix == i ].sum()

    totals = totals.cumsum()
    return arr.cumsum() - totals[ cumsum_ix ]
    # subtract the cumulative previous totalsfor each 

sawtooth1( arr )
# Out[22]: array([ 1,  2,  3,  4,  5, -1, -2, -3, -4,  1,  2,  3,  4, -1, -2, -3])

The time difference on the two is insignificant.

user_8881141 · Answer

Somehow, this looked easier in my head before starting. Luckily, in the meantime, others have already suggested solutions. So, time for comparison of my solution to those of others:

import numpy as np
from timeit import timeit

# original
def f1(signal):
    signal_prev = signal[:-1]
    signal_prev = np.pad(signal_prev,(1,0), mode='constant', constant_values=(0))
    arr = np.array([signal,signal_prev], dtype=int)
    arr = np.transpose(arr)
    
    signal_sum = np.array([], dtype=int)

    for x in arr:
        if np.sign(x[0]*x[1]) > 0:
            signal_sum = np.append(signal_sum, signal_sum[-1] + x[1])
        else:
            signal_sum= np.append(signal_sum, x[0])
    
    arr_sum = np.array([signal, signal_sum])
    return np.transpose(arr_sum)


#Pierre D
def f2(s):
    # make a copy and ensure np.array (in case list was given)
    s = np.array(s).copy()
    idx = np.nonzero(np.diff(s))[0]  # last of each group
    off = np.diff(np.concatenate(([0], np.cumsum(s)[idx])))
    s[idx + 1] -= off
    return np.cumsum(s)


#sai
def f3(signal):
    normal_cumsum = np.cumsum(signal)

    reset_sums = np.roll(np.where(np.diff(np.sign(signal), append=np.sign(signal[-1])) != 0, normal_cumsum, np.zeros_like(signal)), 1)
    ffill_idxs = np.hstack((np.squeeze(np.argwhere(reset_sums != 0)), np.array([len(signal)])))
    
    for start, end in zip(ffill_idxs[:-1], ffill_idxs[1:]):
        reset_sums[start:end] = reset_sums[start]

    return normal_cumsum - reset_sums


#Tis Chris sawtooth 0 
def f4(arr):
    cumsum_ix = np.zeros_like(arr)
    cumsum_ix[ 1: ] = (( arr[ :-1 ] - arr[ 1: ] ) != 0 ).cumsum()
    # cumsum_ix incrememnts for each sign change. 

    result = np.zeros_like( arr )
    for i in range( 0, cumsum_ix[ -1 ] + 1 ):
        # For each cumsum_ix select those items and generate the cumsums.
        result[ cumsum_ix==i ] = arr[ cumsum_ix==i ].cumsum()
    return result
    

#c'est moi
def f5(signal):
    ind = np.where(np.diff(signal, prepend=0).astype(bool))
    signalcount = signal.copy()
    signalcount[ind] = signalcount[ind] + np.diff(ind, prepend=0) * np.sign(signal[ind])
    return signalcount.cumsum()

First, let's check that they all return the same result. Turns out, sai's solution generates sometimes a wrong array:

nxd = 20    #array length    
s = np.random.choice([1, -1], size=nxd)

#Integrity check
print(f1(s).T[1])
#[-1 -2  1 -1  1  2 -1 -2 -3  1  2  3 -1 -2  1 -1  1 -1  1 -1]
print(f2(s))
#[-1 -2  1 -1  1  2 -1 -2 -3  1  2  3 -1 -2  1 -1  1 -1  1 -1]
print(f3(s)) 
#sometimes the counting is incorrect
#[-1 -2  1 -1  1  2  1  ->0<- -1  1  2  3  2  1  1 -1  1 -1  1 -1]
print(f4(s))
#[-1 -2  1 -1  1  2 -1 -2 -3  1  2  3 -1 -2  1 -1  1 -1  1 -1]
print(f5(s))
#[-1 -2  1 -1  1  2 -1 -2 -3  1  2  3 -1 -2  1 -1  1 -1  1 -1]

Now, the timing for shorter arrays:

ntime = 50 #number of test runs
nxd = 2000    #array length

s = np.random.choice([1, -1], size=nxd)

print(timeit(lambda: f1(s), number=ntime))
#1.1028546000000001

print(timeit(lambda: f2(s), number=ntime))
#0.004060500000000022   <--- best timing

print(timeit(lambda: f3(s), number=ntime))
#0.0505989

print(timeit(lambda: f4(s), number=ntime))
#0.4808455999999999

print(timeit(lambda: f5(s), number=ntime))
#0.0046319999999999695

and longer arrays:

ntime = 5 #number of test runs
nxd = int(1e6)    #array length

s = np.random.choice([1, -1], size=nxd)

#print(timeit(lambda: f1(s), number=ntime))
#It took too long to wait for this result

print(timeit(lambda: f2(s), number=ntime))
#0.22104029999999986  <--- again the best timing

print(timeit(lambda: f3(s), number=ntime))
#2.402051

#print(timeit(lambda: f4(s), number=ntime))
#took too long

print(timeit(lambda: f5(s), number=ntime))
#0.2280369000000002

So, the suggestion by Pierre D is your best choice (so far).

How to apply cummulative sum in Numpy in slices with condition to previous value?

4 Answers