Home > Archive > Fortran > September 2005 > pseudoinverse
You are viewing an archived Text-only version of the thread.
To view this thread in it's original format and/or if you want to reply to
this thread please [click here]
|
|
| andrea mordini 2005-09-16, 6:59 pm |
| I'm searching the fortran code to compute the pseudoinverse matrix by
the SVD single value decomposition which i have already computed.
thanks
andrea
| |
| SC Huang 2005-09-16, 6:59 pm |
| OK, I am not sure if this is what you are looking for:
If the svd of matrix is A is:
A = U S V^*
then the pseduoinverse of A is:
pseudo-inv(A) = V S^{-1} U^*
Note that S is diagonal, so S^{-1} is easy to calculate. Also, V^*
denotes the hermitian of V (complex cojugate transpose). If V is real,
then V^* is just the transpose of V. Hope this helps.
On 2005-09-16, andrea mordini <stredostudio@libero.it> wrote:
> I'm searching the fortran code to compute the pseudoinverse matrix by
> the SVD single value decomposition which i have already computed.
> thanks
>
> andrea
| |
| andrea mordini 2005-09-17, 4:03 am |
| What u write is correct but i've tried to compute the SVD
with scilab and with the code provided by numerical recipes
and i obtain different U S and V matrices.
In both the cases A=U S transpose(V) but the calculation of the pinv is
different in my code (wrong) and in scilab (correct)
How is this possible?
andrea
SC Huang wrote:[color=darkred]
> OK, I am not sure if this is what you are looking for:
>
> If the svd of matrix is A is:
> A = U S V^*
> then the pseduoinverse of A is:
> pseudo-inv(A) = V S^{-1} U^*
>
> Note that S is diagonal, so S^{-1} is easy to calculate. Also, V^*
> denotes the hermitian of V (complex cojugate transpose). If V is real,
> then V^* is just the transpose of V. Hope this helps.
>
>
> On 2005-09-16, andrea mordini <stredostudio@libero.it> wrote:
>
| |
| Gordon Sande 2005-09-17, 7:01 pm |
| On 2005-09-17 05:18:12 -0300, andrea mordini <stredostudio@libero.it> said:
> What u write is correct but i've tried to compute the SVD
> with scilab and with the code provided by numerical recipes
> and i obtain different U S and V matrices.
> In both the cases A=U S transpose(V) but the calculation of the pinv is
> different in my code (wrong) and in scilab (correct)
Rude question - How do know which one is right?
Or you may in fact know which one is right and have not yet noticed
your programming error.
> How is this possible?
If the singular values are not well determined then there will be
corresponding ambiguities in the vectors. In particular, if there
are ties then which one is the first one? All that is determined is
the subspace and not its representation. That means that there can
arbitrary orthogonal transformations permitted which means there can
be many right answers, each corresponding to a different choice of
the arbitrary transformation.
[color=darkred]
> andrea
>
>
> SC Huang wrote:
| |
| SC Huang 2005-09-17, 7:01 pm |
| Are the elements in S (diagonal) matrix sorted according to their
alphebatic order in your numerical recipes code? Is it possible that
the numerical recipes code scale S (and either U and V) differently ?
Or is it possible that it returns "V" instead of "transpose(V)" (some
software does this)?
These are just my guessings. I only have experience using matlab and
other fortran codes for svd, but not scilab and numerical receipe's
code. :(
On 2005-09-17, andrea mordini <stredostudio@libero.it> wrote:[color=darkred]
> What u write is correct but i've tried to compute the SVD
> with scilab and with the code provided by numerical recipes
> and i obtain different U S and V matrices.
> In both the cases A=U S transpose(V) but the calculation of the pinv is
> different in my code (wrong) and in scilab (correct)
> How is this possible?
>
> andrea
>
>
> SC Huang wrote:
| |
|
| andrea mordini wrote:
>
> What u write is correct but i've tried to compute the SVD
> with scilab and with the code provided by numerical recipes
> and i obtain different U S and V matrices.
> In both the cases A=U S transpose(V) but the calculation of the pinv is
> different in my code (wrong) and in scilab (correct)
> How is this possible?
With the tools you mentioned - anything goes, and occasionally you might
even get lucky. In the meantime, you need an independent arbitrator.
You'll find it on netlib at,
ieeecss/cascade/mpinv.f
which calculates Moore-Penrose pseudo-inverse of real general matrix.
| |
| andrea mordini 2005-09-19, 3:56 am |
| I'm using now the mpinv.f subroutine
and everything seems to be ok.
thanks to everyone
andrea
| |
| Herman D. Knoble 2005-09-19, 7:58 am |
| Andrea: The Generalized inverse, Moore Penrose inverse, pseudoinverse
for your test matrix can be checked via the code:
http://ftp.cac.psu.edu/pub/ger/fortran/hdk/ginv.f90
(replace the small test driver with your own). This algorithm
was published in the CACM in 1966:-)
Skip Knoble
On Fri, 16 Sep 2005 15:44:20 +0200, andrea mordini <stredostudio@libero.it> wrote:
-|I'm searching the fortran code to compute the pseudoinverse matrix by
-|the SVD single value decomposition which i have already computed.
-|thanks
-|
-|andrea
| |
| Herman D. Knoble 2005-09-19, 7:58 am |
| Andrea, if that is so, then this tests the fact that if the matrix, A,
is square, nonsingular, then the GINV(A) = INVERSE(A). Also if non-square
matrix A has linearly independent columns, then GINV(A) - LeastSquaresInverse(A).
Skip
On Mon, 19 Sep 2005 09:34:34 +0200, andrea mordini <stredostudio@libero.it> wrote:
-|I'm using now the mpinv.f subroutine
-|and everything seems to be ok.
-|thanks to everyone
-|
-|andrea
|
|
|
|
|