CODEWITHSUNDEEP
matlabperformancematrix-multiplicationlarge-data
Heaven
Heaven

Reputation: 118

How to speed up this for-loop code (for large matrix `H_sparse`)?

H_sparse is a large matrix with size 20,000-by-5,000. The matrix-vector product dk = A * Del_H; in the code below is time consuming. How can I speed up this code?

This code is another way to get an equivalent result to the built-in function pinv(H_Sparse) in MATLAB. I think MATLAB uses mex files and bsxfun in pinv, so it's fast.

But in theory the below algorithm is faster:

function PINV_H_Spp = Recur_Pinv_Comp( H_Sparse ) 
         L           = 1;
         H_candidate = H_Sparse(:,L); 
         A           = pinv( H_candidate );
         for L = 1:size( H_Sparse, 2 ) - 1
             L = L + 1; 
             Del_H = H_Sparse(:,L);
             dk    = A * Del_H;
             Ck    =     Del_H - H_candidate * dk; 
             Gk    = pinv( Ck ); 
             A     = A - dk * Gk; 
             A(end+1,:) = Gk; 
             H_candidate(:,end+1) = Del_H;
         end
         PINV_H_Spp = A;

The code can be compared with pinv(H_Sparse), using H_Sparse = rand(20000, 5000) as sample data.

Upvotes: 1

Views: 75

Answers (1)

gnovice
gnovice

Reputation: 125854

A few points of improvement:

  • You can change your loop index to 2:size(H_Sparse, 2) and remove the line L = L + 1;.
  • There's no need to create a separate variable H_candidate, since it's only partitions of H_Sparse. Instead, just index H_sparse accordingly and you'll save on memory.
  • Instead of building your matrix A row-by-row, you can preallocate it and update it using indexing. This usually provides some speed-up.
  • Return A as your output. No need to put it in another variable.

Here's a new version of the code incorporating the above improvements:

function A = Recur_Pinv_Comp(H_Sparse)
  [nRows, nCols] = size(H_Sparse);
  A = [pinv(H_Sparse(:, 1)); zeros(nRows-1, nCols)];
  for L = 2:nCols
    Del_H = H_Sparse(:, L);
    dk = A(1:L-1, :)*Del_H;
    Ck = Del_H - H_Sparse(:, 1:L-1)*dk;
    Gk = pinv(Ck);
    A(1:L-1, :) = A(1:L-1, :) - dk*Gk;
    A(L, :) = Gk;
  end
end

In addition, it looks like your calls to pinv only ever operate on column vectors, so you may be able to replace them with a simple array transpose and scaling by the sum of the squares of the vector (which might speed things up a little more):

Gk  = Ck.'./(Ck.'*Ck);

Upvotes: 1

Related Questions