7 views (last 30 days)

Show older comments

Hi all,

Imagine I have 2 matrices A and B and I need the diagonal elements of (A' * B). If say A and B are both m*n matrices which contains column vectors as follows:

A = [a1, ..., an];

B = [b1, ..., bn];

Then the diagonal elements of (A' * B) can be calculated by the vector products of the paired vectors a1'*b1, a2'*b2, ..., an' * bn.

Now if I perform economy size SVD to A and B:

[ua, sa, va] = svd(A, 'econ');

[ub, sb, vb] = svd(B, 'econ');

Then my question is: can I calculate only the diagonal elements of (A' * B) with SVD components [ua, sa, va] and [ub, sb, vb]? Notice only the diagonal elements, none of the other elements can be computed.

My understanding is no, I can't. Because if recovering (A' * B) by the SVD components, say

A' * B = (ua1 * sa1 * va1')' * (ub1 * sb1 * vb1') + ... + (uan * san * van')' * (ubn * sbn * vbn')

then each component product

(uan * san * van')' * (ubn * sbn * vbn')

is a full m*n matrix, which only recovers part of the whole diagonal line. Summing all n full m*n matrices will recover the matrix (A' * B), but this is not cheaper than computing (A' * B), because it involves computation of all elements of (A' * B). So if SVD is involved, there is no efficient way of getting diagonal elements of (A' * B) by only computing the diagonal elements.

I don't know if I'm right, and please excuse me if it's not so clear.

Many thanks!

Anton Semechko
on 28 Jun 2018

i-th diagonal of A'*B equals dot product of the i-th columns of A and B. Here is a demo showing how to recover only the diagonals of A'*B using SVDs of A and B:

% Generate two M-by-N matrices

M=100;

N=100;

A=randn(M,N);

B=randn(M,N);

% Get SVDs

[Ua,Da,Va]=svd(A,0);

[Ub,Db,Vb]=svd(B,0);

% Pre-multiply D and V'

Va=Da*Va';

Vb=Db*Vb';

% Get diagonal entries of C=A'*B; i-th diagonal element = dot product between i-th columns of A and B

C_diag=zeros(N,1);

for n=1:N

A_n=Ua*Va(:,n); % n-th column of A <--> n-th row of transpose(A)

B_n=Ub*Vb(:,n); % n-th column of B

C_diag(n)=A_n'*B_n;

end

% Verify computed solution

C_diag_ref=diag(A'*B);

E_max=max(abs(C_diag-C_diag_ref));

fprintf('Maximum absolute error: %.10E\n',E_max)

When using this approach, over 100 runs, I got maximum and mean errors of 1.59872115546023e-13 and 1.02029495963052e-13, respectively.

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

Start Hunting!