Finding the number of rows to the next row containing a 1

Hi

I have a column vector of 1s and 0s and I want to find find the number of rows to the next row containing a 1. For example:

A = [0 0 0 1 0 1 1 1 0 0 0 0 0 1]';

I would like the code to return

B = [3 2 1 0 1 0 0 0 5 4 3 2 1 0]';

Is there a vectorized way that this can be done?

Thanks in advance

[KSSV]

But all you have is a column vector...not a matrix.

[andyc100]

Sorry, maybe I didn't ask the question correctly. I'm essentially looking for the simpliest way to do this without looping.

Just as additional information. My column vector is usually in the range of 30,000 rows. Not sure if that would affect the efficiency of some results.

Thanks again.

[Stephen23]

A well-written loop would probably be the simplest and most efficient solution.

[Rik]

@Stephen, Probably, although my current solution beats the current solution by Bruno by a factor of about 10. I don't know what speedup is still possible with either solution.

[Bruno Luong]

"I don't know what speedup is still possible with either solution."

Humm I can not let this challenge passing me.

A=rand(30000,1)>0.7;
ntries = 1000

tic
for i=1:ntries
    % Rik's implementation
    B=A;
    pad=B(end)~=1;
    if pad,B(end+1)=1;end %this method requires the last position to be a 1
    B=flipud(B);
    C=zeros(size(B));
    C(B==1)=[0;diff(find(B))];
    C=ones(size(B))-C;
    out=cumsum(C)-1;
    out=flipud(out);
    if pad,out(end)=[];end
end
trik=toc;

tic
for i=1:ntries
    % Bruno's implementation
    i1=(size(A,1)+1)-flip(find(A));
    B=ones(size(A));
    B(i1)=1-diff([0;i1]);
    B=flip(cumsum(B));
end
tbruno=toc;

fprintf('t Rik   = %f [s]\n', trik);
fprintf('t Bruno = %f [s]\n', tbruno);
fprintf('ratio   = %f\n', trik/tbruno);

Results

t Rik   = 0.725155 [s]
t Bruno = 0.430535 [s]
ratio   = 1.684309 % I BEAT YOU ;-)

[Rik]

I should challenge you more often ;)

[andyc100]

Amazing Rik and Bruno, much appreciate the help.

[Stephen23]

Vectorization is not always the fastest, nor always the clearest code:

A = double(rand(30000,1)>0.7);
%A = [0;0;0;1;0;1;1;1;0;0;0;0;0;1];
N = 1000;

tic
for ii=1:N
	% simple loop
	B0 = A;
	cnt = 0;
	for k = numel(A):-1:1
		cnt = (cnt+1)*(1-A(k));
		B0(k) = cnt;
	end
end
toc

tic
for ii=1:N
	% Rik's implementation
	B=A;
	pad=B(end)~=1;
	if pad,B(end+1)=1;end %this method requires the last position to be a 1
	B=flipud(B);
	C=zeros(size(B));
	C(B==1)=[0;diff(find(B))];
	C=ones(size(B))-C;
	B1=cumsum(C)-1;
	B1=flipud(B1);
	if pad,B1(end)=[];end
end
toc

tic
for ii=1:N
	% Bruno's implementation
	i1=(size(A,1)+1)-flipud(find(A));
	B2=ones(size(A));
	B2(i1)=1-diff([0;i1]);
	B2=flipud(cumsum(B2));
end
toc

isequal(B0,B1,B2)

Giving:

Elapsed time is 0.379344 seconds. % simple loop
Elapsed time is 1.111552 seconds. % Rik
Elapsed time is 0.554412 seconds. % Bruno Luong
ans =
     1

[Bruno Luong]

Very convincing Stephen. In my laptop it's even more obvious

t Rik      = 0.680010 [s]
t Bruno    = 0.378046 [s]
t Stephen  = 0.107799 [s]

In your for-loop code I would cast B initialization in double (in case A is logical)

B = double(A);

or

B = zeros(size(A));

[Bruno Luong]

As I'm slightly surprised by the performance of the for-loop (I would expext it's good but not THAT good), I then try to see how far it from a MEX implementation. And I'm stunned, it's almost as fast (even faster for smaller array).

t Rik      = 0.651531 [s]
t Bruno    = 0.379362 [s]
t Stephen  = 0.104442 [s]
t MEX      = 0.073168 [s]

The code (benchmark + Cmex) is in the attacheh file for those who wants to play with.

I must congratulat TMW for improving the for-loop performance over many years.

Answer 1

[Rik]

It took some time, but here is a solution that should also work for large matrices.

clc,clear
format compact
A = [0 0 0 1 0 1 1 1 0 0 0 0 0 1]';
%    3 2 1 0 1 0 0 0 5 4 3 2 1 0

B=A;
pad=B(end)~=1;
if pad,B(end+1)=1;end %this method requires the last position to be a 1
B=flipud(B);
C=zeros(size(B));
C(B==1)=[0;diff(find(B))];
C=ones(size(B))-C;
out=cumsum(C)-1;
out=flipud(out);
if pad,out(end)=[];end

%only for display:
[A out]

[andyc100]

Thank you. I still don't quite understand the steps (need to think about it a bit more), but it works really well and fast.

[Bruno Luong]

When I run this code with

A = [0;0;0]

I get

out=[3;2;1]

Does it meet the description "find the number of rows to the next row containing a 1".

To me no since there is no 1 in A it should return [0;0;0].

My code also have this flaw.

[Rik]

If that is indeed a problem that should be easy to fix. I chose to write a comment instead of changing that behavior, but maybe I should have made it more explicit. Thank you for drawing more attention to that point.

[andyc100]

This is not a concern for me as there will always be 1s in my implementation. Can always just do a check for not all zeros in the vector before running the code I guess.

Answer 2

[Binbin Qi]

A = [0 0 0 1 0 1 1 1 0 0 0 0 0 1]';
C = find(A);
D = (1:length(A)) - C;
D(D>0) = D(D>0) + inf';
min(abs(D))'

ans =

     3
     2
     1
     0
     1
     0
     0
     0
     5
     4
     3
     2
     1
     0

[andyc100]

Thanks. But I get an error on line:

D = (1:length(A)) - C;

The error is:

Error using  - 
Matrix dimensions must agree.

I forgot to mention that the my column vectors are usually around 30,000 rows. This line looks like it's attempting to create a nxn matrix? Is that correct? If so, will it be slower than a looping method?

[Binbin Qi]

you must use higher version. And this solution need more memory if the vector is large.

[madhan ravi]

 D = bsxfun(@minus, 1:length(A) , C); % for older versions

[Binbin Qi]

if data is larger, you can use the code following

A = [0 0 0 1 0 1 1 1 0 0 0 0 0 1]';
A = [A;1]; % 
B = cell2mat(cellfun(@(x)(length(x)-1:-1:0)',...
    mat2cell(A, diff([0;find(A==1)])),...
    'UniformOutput',false));
B(end) = []

[andyc100]

Thanks so much for this. Always like one liner answers, but have chosen Rik's solution as it works a lot faster for the range of rows I'm working with.

[Andrei Bobrov]

+1

Answer 3

[Bruno Luong]

As much as I love vectorization, this problem is a typical case where the for-loop method is easier, faster, more readable.

This code is ready for 2D array, it works along the first dimension independently.

A=rand(30000,1000)>0.7;

tic
Al=logical(A);
B=zeros(size(A));
b=B(1,:);
for k=size(B,1):-1:1
    b = b+1;
    b(Al(k,:))=0;
    B(k,:)=b;
end
toc

[andyc100]

Thank you Burno. This actually works very fast too and I like that it works for multiple columns. I might end up using it if I end up needing multiple columns.

Answer 4

[Bruno Luong]

Now I just discover CUMSUM has direction option, this is based on Rik's cumsum method, but avoid the double flipping.

B = ones(size(A));
i1 = find(A);
B(i1) = 1-diff([i1;size(A,1)+1]);
B = cumsum(B,'reverse');

[Rik]

Cool, I didn't remember that was an option. Turns out that is already an option as far back as R2015a. The release notes no longer allow you to look back further than R2015a, even if you try modifying the url. All I know is that it isn't possible in R2011a.