Solving Generalized eigenvectors and eigen values to get a common basis

Question

I am looking for solving a generalized eigenvectors and eigen value problem in Matlab. For this, I have tested 2 methods.

if Generalized problem is formulated as :

Then, we could multiply by B^(-1) on each side, such as :

So, from a theorical point of view, it is the simple and classical eigen values problem.

Finally, in Matlab, I simply did, with A=FISH_sp and B=FISH_xc :

[Phi, Lambda] = eig(inv(FISH_xc)*FISH_sp);

But results are not correct when I make after a simple Fisher synthesis (constraints are too bad and also making appear nan values. I don't know why I don't get the same results than the second ones below.

the second method comes from the following paper.

To summarize, the algorithm used is described on page 7. I have followed all the steps of this algorithm and it seems to give better results when I make the Fisher synthesis.

Here the interested part (sorry, I think Latex is not available on stakoverflow) :

Here my little Matlab script for this method :

% Diagonalize A = FISH_sp and B = Fish_xc
[V1,D1] = eig(FISH_sp);
[V2,D2] = eig(FISH_xc);

% Applying each  step of algorithm 1 on page 7
phiB_bar = V2*(D2.^(0.5)+1e-10*eye(7))^(-1);
barA = inv(phiB_bar)*FISH_sp*phiB_bar;
[phiA, vA] = eig(barA);
Phi = phiB_bar*phiA;

So at the end, I find phi eigenvectors matrix (phi) and lambda diagonal matrix (D1).

Now, I would like to do the link between this generalized problem and the eventual common eigenvectors between A and B matrices (respectively Fish_sp and Fish_xc). Is there a way to perform this?

Indeed, what I have done up to now is to to find a parallel relation between A*Phi and B*Phi, linked by Lambda diagonal matrix. Maybe, we could arrange this relation such that :

A*Phi'=Phi'*Lambda_A'

and

B*Phi'=Phi'*Lambda_B'

From a numerical point of view, why don't I get the same results between the method in 1) and the method in 3) ? I mean about the Phi eigen vectors matrix and the Lambda diagonal matrix.

However, this is the same formulation.

EDIT :

I have wrong results if I want to say that phi diagonalizes both A=FISH_sp and B=FISH_xc matrices.

Indeed, by doing :

% Marginalizing over uncommon parameters between the two matrices
COV_GCsp_first = inv(FISH_GCsp);
COV_XC_first = inv(FISH_XC);
COV_GCsp = COV_GCsp_first(1:N,1:N);
COV_XC = COV_XC_first(1:N,1:N);
% Invert to get Fisher matrix
FISH_sp = inv(COV_GCsp);
FISH_xc = inv(COV_XC);
% Diagonalize
[V1,D1] = eig(FISH_sp);
[V2,D2] = eig(FISH_xc);

% Build phi matrix
% V2 corresponds to eigen vectors of FISH_xc
phiB_bar = V2*diag(diag(D2.^(-0.5)));
% DEBUG : check identity matrix => OK, Identity matrix found !
id = (phiB_bar')*FISH_xc*phiB_bar
% phi matrix
barA = (phiB_bar')*FISH_sp*phiB_bar
[phiA, vA] = eig(barA);
phi = phiB_bar*phiA;

% Check eigen values : OK, columns of eigenvalues found !
FISH_sp*V1./V1
% Check eigen values : OK, columns of eigenvalues found !
FISH_xc*V2./V2

% Check if phi diagolize FISH_sp : NOT OK, not identical eigenvalues 
FISH_sp*phi./phi
% Check if phi diagolize FISH_sp : NOT OK, not identical eigenvalues 
FISH_xc*phi./phi

So, I don't find that matrix of eigenvectors Phi diagonalizes A and B since the eigenvalues expected are not columns of identical values.

By the way, I find the eigenvalues D1 and D2 coming from :

[V1,D1] = eig(FISH_sp);
[V2,D2] = eig(FISH_xc);

% Check eigen values : OK, columns of eigenvalues D1 found !
FISH_sp*V1./V1
% Check eigen values : OK, columns of eigenvalues D2 found !
FISH_xc*V2./V2

How could I fix this wrong result (I am talking about the ratios :

FISH_sp*phi./phi
FISH_xc*phi./phi

which don't give same values for a given column of FISH_sp and FISH_xc) ) ?? In the paper, they say that phi diagonalizes A=FISH_sp and B=FISH_xc but I can't reproduce it.

If someone could see where is my error ...

I only quickly read through half of the question, but I can already see an error in your definition phiB_bar = V2*(D2.^(0.5)+1e-10*eye(7)^(-1)); in your second method: from the formula from the paper that should be phiB_bar = V2*inv(D2.^(0.5)+1e-10*eye(7));, shouldn't it? — Matteo V
@MatteoV . Thanks for this suggestion, you are right : I have replace by phiB_bar = V2*(D2.^(0.5)+1e-10*eye(7))^(-1); But the phi doesn't still appear as a eigenvector of both A=FISH_sp and B=FISH_xp. It is frustrating... — youpilat13
Just did a quick test with two 5x5 random matrices, and the method 1 seems to be correct (norm(A*phi - B*phi*lambda) = 1.933967216188714e-14 in my case), while for method 2 it returns norm(A*phi - B*phi*lambda) = 2.228744851244737. Maybe there are some underlying assumptions in the paper (disclaimer: I haven't read it)? — Matteo V
[V, D] = eig(A, B) should solve the Generalized Eigenvalue Problem A*V = B*V*D. Is this not working for your problem? For questions about the paper math.stackexchange.com may be the better place to ask. — pH 74
@pH74 yes, I can use directly the computation with [V, D] = eig(A, B) but I would like to do the link with a common eigenevctors basis between A and B. One suggested me to look a the problem of joint diagonalization, and it has two variants: orthogonal, in which the basis vectors are orthonormal, and non-orthogonal, which is harder to solve, but which may be more appropriate to your application. You can see this advise on : scicomp.stackexchange.com/questions/36670/… . — youpilat13

Luk Luk · Accepted Answer · 2021-02-02T14:42:45

You defined barA = inv(phiB_bar)*FISH_sp*phiB_bar. From Eq. (39) in the manuscript it looks like it should be barA = transpose(phiB_bar)*FISH_sp*phiB_bar instead.

Also, your mehtod 1 fails when B is singular (an inverse does not exist). MATLAB's eig(A,B) however should handle also singular Bs if memory serves me well.

Solving Generalized eigenvectors and eigen values to get a common basis

1 Answers