cumsum with upper & lower limits?

Question

I'd like to find a vectorized way to calculate the cumulative sums of a vector, but with upper and lower limits.

In my case, the input only contains 1's and -1's. You can use this assumption in your answer. Of course, a more general solution is also welcome.

For example:

x     = [1 1 1 1 -1 -1 -1 -1 -1 -1];
upper = 3;
lower = 0;
s     = cumsum(x)                    %// Ordinary cumsum.
s =
     1     2     3     4     3     2     1     0    -1    -2

y     = cumsumlim(x, upper, lower)   %// Cumsum with limits.
y =
     1     2     3     3     2     1     0     0     0     0
                 ^                       ^
                 |                       |
            upper limit             lower limit

When the cumulative sum reaches the upper limit (at the 3rd element), it won't increase anymore. Likewise, when the cumulative sum reaches the lower limit (at the 7th element), it won't decrease anymore. A for-loop version would be like this:

function y = cumsumlim(x, upper, lower)

y = zeros(size(x));
y(1) = x(1);

for i = 2 : numel(x)
    y(i) = y(i-1) + x(i);
    y(i) = min(y(i), upper);
    y(i) = max(y(i), lower);
end

end

Do you have any ideas?

Answer 1

This is a somewhat hackish solution, but perhaps worth mentioning.

You can do the sum using a signed integer data type, and exploit the inherent limits of that data type. For this to work, the input needs to be converted to that integer type and multiplied by the appropiate factor, and an initial offset needs to be applied. The factor and offset are chosen as a function of lower and upper. After cumsum, the multiplication and offset are undone to obtain the desired result.

In your example, data type int8 suffices; and the required factor and offset are 85 and -128 respectively:

x = [1 1 1 1 -1 -1 -1 -1 -1 -1];
result = cumsum([-128 int8(x)*85]); %// integer sum, with factor and initial offset
result = (double(result(2:end))+128)/85; %// undo factor and offset

which gives

result =
     1     2     3     3     2     1     0     0     0     0

Answer 2

I won't provide you with a magic vectorized way to do this, but I'll provide you with some data that probably will help you get on with your work.

Your cumsumlim function is very fast!

tic
for ii = 1:100
    y = cumsumlim(x,3,0);
end
t = toc;
disp(['Length of vector: ' num2str(numel(x))])
disp(['Total time for one execution: ' num2str(t*10), ' ms.'])
Length of vector: 65000
Total time for one execution: 1.7965 ms.

I really doubt this is your bottleneck. Have you tried profiling the code?