Pourquoi ne puis-je pas obtenir un SVD valide de X via la décomposition de valeurs propres de XX 'et X'X?

J'essaie de faire SVD à la main:

m<-matrix(c(1,0,1,2,1,1,1,0,0),byrow=TRUE,nrow=3)

U=eigen(m%*%t(m))$vector
V=eigen(t(m)%*%m)$vector
D=sqrt(diag(eigen(m%*%t(m))$values))

U1=svd(m)$u
V1=svd(m)$v
D1=diag(svd(m)$d)

U1%*%D1%*%t(V1)
U%*%D%*%t(V)

Mais la dernière ligne ne revient pas m. Pourquoi? Cela semble avoir quelque chose à voir avec les signes de ces vecteurs propres ... Ou ai-je mal compris la procédure?

Analyse du problème

La SVD d'une matrice n'est jamais unique. Soit la matrice dimensions et que sa SVD soitn × $A$$n\times k$

$$A = U D V^\prime$$

pour un la matrice avec des colonnes orthonormales, une diagonale la matrice avec des entrées non-négatives, et un la matrice avec des colonnes orthonormales.U p × p D k × p $n\times p$$U$$p\times p$$D$$k\times p$$V$

Maintenant, choisissez arbitrairement n'importe quelle diagonale matrice ayant s sur la diagonale, de sorte que soit l' identité . alorsS ± 1 S 2 = I p × p I $p\times p$$S$$\pm 1$$S^2 = I$$p\times p$$I_p$

$$A = U D V^\prime = U I D I V^\prime = U (S^2) D (S^2) V^\prime = (US) (SDS) (VS)^\prime$$

est également un SVD de parce que montre que a des colonnes orthonormées et un calcul similaire montre que a des colonnes orthonormées. De plus, puisque et sont diagonaux, ils font la navette, d'où montre que toujours des entrées non négatives.( U S ) ′ ( U S ) = S ′ U ′ U S = S ′ I p S = S ′ S = S 2 = I p U S V S S D S D S = D S 2 = D $A$

$$(US)^\prime(US) = S^\prime U^\prime U S = S^\prime I_p S = S^\prime S = S^2 = I_p$$ $US$$VS$$S$$D$ $$S D S = DS^2 = D$$ $D$

La méthode implémentée dans le code pour trouver un SVD trouve un qui diagonalise et, de la même façon, un qui diagonale Il procède au calcul de en fonction des valeurs propres trouvées dans . Le problème est que cela ne garantit pas une correspondance cohérente des colonnes de avec les colonnes de .A A ′ = ( U D V ′ ) ( U D V ′ ) ′ = U D V ′ V D ′ U ′ = U D 2 U ′ V A ′ A = V D 2 V ′ . D D 2 U V$U$

$$AA^\prime = (UDV^\prime)(UDV^\prime)^\prime = UDV^\prime V D^\prime U^\prime = UD^2 U^\prime$$ $V$ $$A^\prime A = VD^2V^\prime.$$ $D$$D^2$$U$$V$

Une solution

Au lieu de cela, après avoir trouvé un tel et un tel , utilisez-les pour calculer$U$$V$

$$U^\prime A V = U^\prime (U D V^\prime) V = (U^\prime U) D (V^\prime V) = D$$

directement et efficacement. Les valeurs diagonales de ce ne sont pas nécessairement positives. $D$ (C'est parce qu'il n'y a rien dans le processus de diagonalisation de ou qui garantira que, puisque ces deux processus ont été effectués séparément.) Rendez-les positifs en choisissant les entrées le long de la diagonale de pour égaler les signes des entrées de , de sorte que ait toutes les valeurs positives. Compensez cela en multipliant à droite par :$A^\prime A$$AA^\prime$$S$$D$$SD$$U$$S$

$$A = U D V^\prime = (US) (SD) V^\prime.$$

C'est un SVD.

Exemple

Soit avec . Un SVD est$n=p=k=1$$A=(-2)$

$$(-2) = (1)(2)(-1)$$

avec , et .$U=(1)$$D=(2)$$V=(-1)$

Si vous diagonalisez$A^\prime A = (4)$$U=(1)$$D=(\sqrt{4})=(2)$$AA^\prime=(4)$$V=(1)$

$$UDV^\prime = (1)(2)(1) = (2) \ne A.$$ $$D=U^\prime A V = (1)^\prime (-2) (1) = (-2).$$ $S=(-1)$$U$$US = (1)(-1)=(-1)$$D$$SD = (-1)(-2)=(2)$ $$A = (-1)(2)(1),$$

Code

Voici le code modifié. Sa sortie confirme

1. La méthode recrée mcorrectement.
2. $U$ et sont vraiment encore orthonormés.$V$
3. Mais le résultat n'est pas le même SVD renvoyé par svd. (Les deux sont également valables.)

m <- matrix(c(1,0,1,2,1,1,1,0,0),byrow=TRUE,nrow=3)

U <- eigen(tcrossprod(m))$vector
V <- eigen(crossprod(m))$vector
D <- diag(zapsmall(diag(t(U) %*% m %*% V)))
s <- diag(sign(diag(D)))  # Find the signs of the eigenvalues
U <- U %*% s              # Adjust the columns of U
D <- s %*% D              # Fix up D.  (D <- abs(D) would be more efficient.)

U1=svd(m)$u
V1=svd(m)$v
D1=diag(svd(m)$d,n,n)

zapsmall(U1 %*% D1 %*% t(V1)) # SVD
zapsmall(U %*% D %*% t(V))    # Hand-rolled SVD
zapsmall(crossprod(U))        # Check that U is orthonormal
zapsmall(tcrossprod(V))       # Check that V' is orthonormal

Comme je l'ai souligné dans un commentaire à la réponse de @ whuber, cette méthode pour calculer le SVD ne fonctionne pas pour chaque matrice . Le problème ne se limite pas aux panneaux.

Le problème est qu'il peut y avoir des valeurs propres répétées, et dans ce cas, la composition propre de et n'est pas unique et tous les choix de et ne peuvent pas être utilisés pour récupérer le facteur diagonal de la SVD. Par exemple, si vous prenez une matrice orthogonale non diagonale (disons ), alors . Parmi tous les choix possibles pour la matrice de vecteur propre de , retournera , donc dans ce cas n'est pas diagonal.A A ' U V A = [ 3 / 5 4 / 5 - 4 / 5 3 / 5 ] A A ' = A ' A = I I U = V = I U ' A V = $A'A$$AA'$$U$$V$$A=\begin{bmatrix}3/5&4/5\\-4/5&3/5\end{bmatrix}$$AA'=A'A=I$$I$eigen$U=V=I$$U'AV=A$

Intuitivement, ceci est une autre manifestation du même problème que @whuber souligne, qu'il doit y avoir une "correspondance" entre les colonnes de et , et le calcul de deux compositions à eigend séparément ne le garantit pas.$U$$V$

Si toutes les valeurs singulières de sont distinctes, alors la composition par eigend est unique (jusqu'à l'échelle / signes) et la méthode fonctionne. Remarque: ce n'est toujours pas une bonne idée de l'utiliser en code de production sur un ordinateur avec une arithmétique à virgule flottante, car lorsque vous les produits et le résultat calculé peut être perturbé par une quantité de l'ordre de , où est la précision de la machine. Si les grandeurs des valeurs singulières diffèrent considérablement (de plus de , grosso modo), cela nuit à la précision numérique des plus petites.A ′ A A A ′ ‖ A ‖ 2$A$$A'A$$AA'$$\|A\|^2u$$u \approx 2\times 10^{-16}$$10^{-8}$

Le calcul de la SVD à partir des deux compositions de type eigend est un excellent exemple d'apprentissage, mais dans les applications réelles, utilisez toujours la svdfonction de R pour calculer la décomposition en valeurs singulières.