Cumsum reversed?

General Usage
1

I need to do cumsum where I start at the end of vector v, and sum backwards.

Is there a separate function for this, or should I just do cumsum(v[end:-1:1])[end:-1:1]?

1 Like
2

Not pretty but easy to remember: reverse(cumsum(reverse(v)))

NB:
And for a faster in-place version of the above:
reverse!(cumsum(reverse!(x)))

3 Likes
3

If you want performance you probably need to use a custom function. This should be close to optimal:

julia> function rev_cumsum!(c,x)
           cumsum = zero(eltype(x))
           @inbounds for i in lastindex(x):-1:firstindex(x)
               cumsum += x[i]
               c[i] = cumsum
           end
           return c
       end
rev_cumsum! (generic function with 1 method)

julia> function rev_cumsum(x)
           c = similar(x)
           return rev_cumsum!(c,x)
       end
rev_cumsum (generic function with 1 method)

julia> x = rand(1:10,1000);

julia> @btime reverse(cumsum(reverse($x)));
  2.516 μs (3 allocations: 23.81 KiB)

julia> @btime rev_cumsum($x);
  570.379 ns (1 allocation: 7.94 KiB)

julia> @btime rev_cumsum!(c,$x) setup=(c=similar(x));
  364.524 ns (0 allocations: 0 bytes)

(edit: added @inbounds, which increases significantly the performance)

2 Likes
4

Maybe worth noting that you can use Iterators.reverse for this, at least for vectors, although you still need to reverse the result. You can also make a view with reversed indices, which will work for cumsum(x; dims) on matrices etc. too. Should be comparable to the hand-written loop:

julia> @btime reverse(cumsum(reverse($x)));
  min 1.812 μs, mean 5.976 μs (3 allocations: 23.81 KiB, mean GC time 22.64%)

julia> @btime rev_cumsum($x);
  min 424.731 ns, mean 1.569 μs (1 allocation: 7.94 KiB, mean GC time 25.95%)

julia> @btime reverse!(cumsum(Iterators.reverse($x)));
  min 736.321 ns, mean 1.882 μs (1 allocation: 7.94 KiB, mean GC time 22.86%)

julia> function rev_cumsum_v(x::AbstractVector)
         y = similar(x)
         cumsum!(view(y, reverse(axes(x,1))), view(x, reverse(axes(x,1))))
         y
       end;

julia> @btime rev_cumsum_v($x);
  min 520.401 ns, mean 1.439 μs (7 allocations: 8.31 KiB, mean GC time 25.81%)
2 Likes
5

There could be a cumsum(f,x) that accepts a function as a first argument, or this is too rare as a use case to be worth the additional complication?

6

I think that would probably be in line with what was done to many of the other statistics function in order to allow the same, yet I don’t think it would help in this case as f is applied to all elements.

1 Like
7

What I had in mind was a mapping, such cumsum(i -> x[end-i+1], x), but that would be something else, of course.