3
$\begingroup$

I'm trying to follow Prof Strang exercise but I have a problem with the signs of the resultant matrix.

He also had a problem with the signs during the exercise, so I have no way to find what is the error.

Can somebody point me what I'm doing wrong?

Thanks in advance,

Diego

These are my results:

> A
     [,1] [,2]
[1,]    4    4
[2,]   -3    3
> SVD_of_A
     [,1] [,2]
[1,]   -4    4
[2,]   -3   -3
> U
     [,1] [,2]
[1,]   -1    0
[2,]    0   -1
> sigma
     [,1]     [,2]
[1,] 5.656854 0.000000
[2,] 0.000000 4.242641
> V
      [,1]       [,2]
[1,] 0.7071068 -0.7071068
[2,] 0.7071068  0.7071068

These are the results of the svd function

> svd(A)
$d
[1] 5.656854 4.242641

$u
     [,1] [,2]
[1,]   -1    0
[2,]    0    1

$v
       [,1]       [,2]
[1,] -0.7071068 -0.7071068
[2,] -0.7071068  0.7071068

This is the code I'm using for that exercise:

rm(list=ls())

#------------------------------------------------------------------
# Find the Singular Value Decomposition of the matrix A
#------------------------------------------------------------------

A = matrix(c(4,-3,4,3),2)

#------------------------------------------------------------------
# We want to obtain:
# SVD = U*sigma*t(V)
#------------------------------------------------------------------

#------------------------------------------------------------------
# First we need to obtain the transpose of A
# t = Given a matrix or data.frame x, t returns the transpose of x.
#------------------------------------------------------------------

At= t(A)

#------------------------------------------------------------------
# Multiply the 2 matrices, At and A to obtain AtA
# %*% = Matrix Multiplication
#------------------------------------------------------------------

AtA=At%*%A

#------------------------------------------------------------------
# Now we need to obtain the spectral decomposition of the matrix AtA
# eigen = Computes eigenvalues and eigenvectors of real 
#         (double, integer, logical) or complex matrices.
#------------------------------------------------------------------

eig1 = eigen(AtA)

#------------------------------------------------------------------
# With the eigenvalues of AtA we can build sigma.
# I did this manually here for the purpose of the procedure.
# However note that sigma is equal to sqrt(AAt)
#------------------------------------------------------------------

sigma=matrix(c(sqrt(eig1$values[1]),0,0,sqrt(eig1$values[2])),2)

#------------------------------------------------------------------
# The eigenvectors of AtA are the values of V
#------------------------------------------------------------------

V=as.matrix(eig1$vectors)

#------------------------------------------------------------------
# Next we need to compute A*At
#------------------------------------------------------------------

AAt=A%*%At

#------------------------------------------------------------------
# Now we need to obtain the spectral decomposition of the matrix AAt
#------------------------------------------------------------------

eig2=eigen(AAt)

#------------------------------------------------------------------
# U is equal to the eigen vectors of AAt
#------------------------------------------------------------------

U=as.matrix(eig2$vectors)

#------------------------------------------------------------------
# Finally the singular value decomposition of A is equal to U*sigma*V
#------------------------------------------------------------------

SVD_of_A=U%*%sigma%*%V

#------------------------------------------------------------------
# Results
#------------------------------------------------------------------

A
SVD_of_A
U
sigma
V

PS: I've also tried V transpose, but the results are still wrong...

$\endgroup$
3
  • 5
    $\begingroup$ You're not doing anything wrong, nor is the software. The SVD is not unique. If we perturb the column vectors of both of the basis matrices by the same sign matrix, we obtain another valid SVD for the same matrix. This is also true for eigendecompositions of symmetric matrices. I didn't run your code, but what almost undoubtedly happened was that each of your eigen calls ended up choosing a different sign pattern for the basis vectors. There's no algorithmic error here. $\endgroup$ – cardinal Feb 6 '13 at 14:45
  • 2
    $\begingroup$ For an extended illustration of this phenomenon (and some ways to deal with it) please see stats.stackexchange.com/questions/34396/…. In addition to what @Cardinal said, there are also the rare cases where the dimension of an eigenspace is greater than $1$: the software is then free to pick any arbitrary orthonormal basis of such an eigenspace. It would be strange for two software solutions ever to agree in such instances, unless they used exactly the same underlying code. $\endgroup$ – whuber Feb 6 '13 at 17:43
  • $\begingroup$ Thank you very much cardinal and whuber for the explanations! I'll take a look at the link. Just wondering... Why or better How the SVD function in R is able to produce always the same result? $\endgroup$ – Diego Feb 7 '13 at 1:01
2
$\begingroup$

I did the exact same thing with Matlab/Octave and got the right result.

A = [4 4;-3 3]

At= transpose(A)
AtA = At * A
AAt = A * At

[vec1, val1] = eig(AtA)
V = vec1

[vec2, val2] = eig(AAt)
U = vec2

sigma=sqrt(val1)

SVD_of_A = U * sigma * transpose(V)

Result:

SVD_of_A =

   4.0000   4.0000
  -3.0000   3.0000

The error is that $eigen(AAt)$ should produce eigenvector equals to have the right result to:

      [,1] [,2]
[1,]    1    0
[2,]    0    1

See comment of Cardinal

$\endgroup$
1
  • 4
    $\begingroup$ No, $R$ did not do anything "wrong" in this instance. This comes down to understanding that the SVD is not unique. The vectors are invariant to common sign flips, and this is precisely what happened here. $\endgroup$ – cardinal Feb 6 '13 at 14:41

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.