Clear Filters
Clear Filters

How to get around sparse row deletion for least squares calculation

4 views (last 30 days)
I have an alrogithm that repeats a reduced basis calculation for a sparse block digaonal matrix Psi. Here Psi is composed of d number of equally sized matrices . I perform a least squares calculation for this matrix with a vector myVec to get the coefficient vector c. I then perforom an algorithm to delete a row from each of in the matrix Psi (as well as from entries in myvec). The problem I'm running into is that this is slow as I'm essentially asking for a new sparse block diagonal matrix of for each iteration, to the point where doing the least squares calculation is faster using a full matrix for Psi. Is there potentially some clever trick I can employ that allows me to delete rows from the sparse matrix without actually changing the size?

Answers (2)

Matt J
Matt J on 6 Mar 2023
Edited: Matt J on 6 Mar 2023
No, probably not. You might try pagemldivide instead, thereby avoiding sparse matrices altogether.
  1 Comment
Matt J
Matt J on 6 Mar 2023
Edited: Matt J on 6 Mar 2023
With the speed-up from pagemldivide, you might be able to offset the cost of rebuilding the matrices every time:
Psi=repmat({rand(30)},1,500); blocks=cellfun(@sparse,Psi,'uni',0);
A=cat(3,Psi{:}); b=A(:,1,:); %paged form (full)
S=blkdiag(blocks{:}); %block diagonal form (sparse)
timeit( @()S\b(:) )
ans = 0.0087
timeit( @() pagemldivide(A,b) )
ans = 0.0015

Sign in to comment.


Bruno Luong
Bruno Luong on 6 Mar 2023
Edited: Bruno Luong on 6 Mar 2023
I'm nit sure why you formalize as block sparse, since it is like solving d independent linear systems of the same size, and can be performed by pagemrdivide as Matt has pointed out.
Now back too your question, you could probably reformulate the row-deletion to a weigted lsq system:
m = 10;
n = 3;
p = 1;
format long
% Full solution
A = rand(m,n);
b = rand(m,p);
x = A\b
x = 3×1
0.092568593269531 0.674057731207846 0.304093420391084
% solution after removed 10th row
xd = A(1:m-1,:)\b(1:m-1,:)
xd = 3×1
-0.215079045247037 0.946905459896963 0.286462609971429
% weigthed least-square solution w
w = ones(m,1);
w (m) = 0;
xw = (w.*A)\(w.*b)
xw = 3×1
-0.215079045247037 0.946905459896963 0.286462609971429
  3 Comments
Bruno Luong
Bruno Luong on 6 Mar 2023
Edited: Bruno Luong on 6 Mar 2023
OP asks to keep the original SPARSE matrix, so I propose a solution for him to avoid : "I'm essentially asking for a new sparse block diagonal matrix of for each iteration".
Furthermore he could use sparse algorithms that requires matrix product vector so instead of doing (w.*A)*x, it only needs w.*(A*x).

Sign in to comment.

Categories

Find more on Operating on Diagonal Matrices in Help Center and File Exchange

Products


Release

R2022a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!