How to speed up matrix update zeros by using previous value?

Hello,

I am looking for the way to update my matrix zero values by using previous value. For example if c matrix is:

0 1 5 2 3 9 0
1 0 5 0 0 7 0
0 0 0 1 7 8 9
7 2 1 0 0 0 7

The result should look like this:

1 1 5 2 3 9 9
1 1 5 2 3 7 9
1 1 5 1 7 8 9
7 2 1 1 7 8 7

Now I use the following code:

c(c == 0) = NaN;
c = fillmissing(c,'previous');
c = fillmissing(c,'nearest',1);

However I use quite big matrixes like 10000000x150 and it takes very long to process, maybe somebody knows of a way to speed up this?

Answer 1

[Jan]

Start with:

c = fillmissing(c, 'previous', 'EndValues', 'nearest');

Another idea is to process the matrix column-wise:

c = [0 1 5 2 3 9 0
     1 0 5 0 0 7 0
     0 0 0 1 7 8 9
     7 2 1 0 0 0 7];

for col = 1:size(c, 2)
    x = c(:, col);         % Get current column
    f = (x ~= 0);          % Find non-zero elements
    m = cumsum(f);         % Cumulative sum of logical indices
    ind1 = find(f, 1);     % Fill initial zeros:
    if ~isempty(ind1)
      m(1:ind1) = 1;
    end
    xf = x(f);             % Vector of non-zero elements
    c(:, col) = xf(m);     % Replace column
end

This avoids the useless replacing of 0 to NaN and the creation of the huge temporary matrices.

This code fails, if the column contains zeros only.

A C-mex file could avoid the creation of f, m and xf. Do you have a C-compiler?

[Mantas Vaitonis]

Thank You for the help. Unfortunately I do not have C-compiler. The second option with for loop seemed to be the most effective.

[Bruno Luong]

If I'm not mistaken, the biggest problem of Jan's code is still in this line

c(:, col) = xf(m);     % Replace column

The entire matrix is still duplicated before the column is replaced (copy on write), so it will have 150 times memory moving.

Copying f, m and xf is not a big deal compared to that.

[Jan]

@Bruno: Are you sure? After the first c(:, 1) = xf(m), this function has its own copy of c . Is there really a need to further deep copies?

I cannot run it currently, but it should be easy to test using:

format debug
c = [0 1 5 2 3 9 0
     1 0 5 0 0 7 0
     0 0 0 1 7 8 9
     7 2 1 0 0 0 7];
c
for col = 1:size(c, 2)
    x = c(:, col);         % Get current column
    f = (x ~= 0);          % Find non-zero elements
    m = cumsum(f);         % Cumulative sum of logical indices
    ind1 = find(f, 1);     % Fill initial zeros:
    if ~isempty(ind1)
      m(1:ind1) = 1;
    end
    xf = x(f);             % Vector of non-zero elements
    c(:, col) = xf(m);     % Replace column
    c
end

If this does really performs a deep copy of c in each iteration, this would be a waste of time. Then:

function R = FillMyGaps(c)
ncol = size(c, 2);
Result = cell(1, ncol)
for col = 1:ncol
    x = c(:, col);         % Get current column
    f = (x ~= 0);          % Find non-zero elements
    m = cumsum(f);           % Cumulative sum of logical indices
    ind1 = find(f, 1);     % Fill initial zeros:
    if ~isempty(ind1)
      m(1:ind1) = 1;
    end
    xf = x(f);             % Vector of non-zero elements
    Result{col} = xf(m);   % Replace column
end
R = cat(2, Result{:});
% Or maybe better: <https://www.mathworks.com/matlabcentral/fileexchange/28916-cell2vec>
% R = reshape(Cell2Vec(Result), [], ncol);

[Bruno Luong]

Sorry Jan, there is indeed of course no deep copy. My bad.

[Jan]

Your interjections are always welcome, because they provide deeper insights very frequently.

I'm still disappointed, that xf=x(f); y=xf(m) cannot be abbreviated. Such a "cumulated indexing" should work in one step.

Nevertheless, in a C-mex, the creation of the vectors x, f, m and xf can be avoided. @Mantas Vaitonis: If this is time-critical, install a C-compiler and ask for the posting of the small function.

[Bruno Luong]

There are shorter code using INTERP1 but both are slighly slower than your code on my benchmark. Also both can work on multicolumns so one can tune up the chunk size depending on the amount of memory available.

for col = 1:size(c, 2)
    x = c(:, col);
    f = (x ~= 0); 
    x(~f) = interp1(find(f),x(f),find(~f),'previous','extrap');
    c(:, col) = x;  
end

Or

for col = 1:size(c, 2)
    x = c(:, col);
    f = (x ~= 0);
    c(:, col) =  interp1(find(f),x(f),1:length(f),'previous','extrap');
end