[R] spectral decomposition for near-singular pd matrices

[Moshe Olshansky]

Your problem seems to be pretty hard (for larger dimensions)!

When speaking about the A1,A2,...,An sequence I did not mean representing A as a product of A1*A2*...*An. What I mean is that A1,...,An is a sequence of matrices which starts at a "good" matrix, i.e. it's eigenvalues are well separated, then A2 is not too different from A1 and having already found the eigensystem (or a part of it) for A1 you can use them to find that system for A2 (without doing it from scratch), A3 is not too much different from A2 and so having the eigensystem for A2 you get ne for A3, etc. until you get to the last matrix which is A.

We had a matrix arising from elasticity. When the elasticity coefficient mu was small the matrix was good but when mu approached 0.5, half of the eigenvalues approached 0 (for mu = 0.5 half of the eigenvalues are really 0).
Starting from something like mu = 0.4 we were able to track several selected (lowest) eigenvalues up to something like mu = 0.4999, while no other program (at that time) could find the eigensystem of such matrix even for much less problematic values of mu.


--- On Thu, 17/7/08, Prasenjit Kapat <kapatp at gmail.com> wrote:

> From: Prasenjit Kapat <kapatp at gmail.com>
> Subject: Re: [R] spectral decomposition for near-singular pd matrices
> To: r-help at stat.math.ethz.ch
> Received: Thursday, 17 July, 2008, 10:56 AM
> Moshe Olshansky <m_olshansky <at> yahoo.com>
> writes:
> 
> > How large is your matrix?
> 
> Right now I am looking at sizes between 30x30 to 150x150,
> though it will 
> increase in future.
> 
> > Are the very small eigenvalues well separated?
> > 
> > If your matrix is not very small and the lower
> eigenvalues are clustered, 
> this may be a really hard problem!
> 
> Here are the deciles of the eigenvalues (using eigen())
> from a typical random 
> generation (30x30):
>           Min          10th         20th         30th      
>   40th
> -1.132398e-16 -6.132983e-17 3.262002e-18 1.972702e-17
> 8.429709e-17
>          50th         60th          70th         80th     
> 90th     Max
>  2.065645e-13  1.624553e-09 4.730935e-06   0.00443353
> 0.9171549 16.5156
> 
> I am guessing they are *not* "well-separated." 
> 
> > You may need a special purpose algorithm and/or higher
> precision arithmetic.
> > If your matrix is A and there exists a sequence of
> matrices A1,A2,...An = A 
> such that A1 is "good", A2 is a bit
> > worse (and is not very different from A1), etc., you
> may be able to compute 
> the eigensystem for A1 and then
> > track it down to A2, A3,...,An. I have worked on such
> problems some 15 years 
> ago but I believe that by now
> > there should be packages doing this (I do not know
> whether they exist in R).
> 
> I will have to think on the possibility of such a
> product-decomposition. But 
> even if it were possible, it would be an infinite product
> (just thinking out 
> loud).
> 
> Thanks again,
> PK
> 
> PS: Kindly CC me when replying.
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained,
> reproducible code.