来自 numpy.linalg.eig 特征向量的任意符号导致自定义 SVD 实现的输出不正确

Question

对于一个类（阅读：不能只使用

np.linalg.svd

）我正在尝试编写一个函数来生成任意形状矩阵的奇异值分解（SVD）；它应该将一些输入矩阵 [A] 分解为 3 个矩阵 [U]、[S] 和 [V]，使得 [A] = [U][S][V]^T 。好的，这是我的尝试，分别从[A][A]T和[A]^T[A]的特征向量得到[U]^和[V]：

import numpy as np
from math import sqrt

def SVD(tIn):
    
    tLeftEigenvalues_unsorted, tLeftColumnEigenvectors_unsorted = np.linalg.eig(tIn @ tIn.T)
    tLeftRowEigenvectors_unsorted = tLeftColumnEigenvectors_unsorted.T
    tDescendingIndices = np.argsort(tLeftEigenvalues_unsorted)[::-1]
    tLeftEigenvalues = np.array(tLeftEigenvalues_unsorted)[tDescendingIndices]
    tLeftRowEigenvectors = np.array(tLeftRowEigenvectors_unsorted)[tDescendingIndices]

    tRightEigenvalues_unsorted, tRightColumnEigenvectors_unsorted = np.linalg.eig(tIn.T @ tIn)
    tRightRowEigenvectors_unsorted = tRightColumnEigenvectors_unsorted.T
    tDescendingIndices = np.argsort(tRightEigenvalues_unsorted)[::-1]
    tRightEigenvalues = np.array(tRightEigenvalues_unsorted)[tDescendingIndices]
    tRightRowEigenvectors = np.array(tRightRowEigenvectors_unsorted)[tDescendingIndices]
      
    tSingularValues = [sqrt(i) for i in tRightEigenvalues]
    
    tS = np.zeros(tIn.shape)
    np.fill_diagonal(tS, tSingularValues)
    
    tU = tLeftRowEigenvectors.T
    tVt = tRightRowEigenvectors
    
    return tU, tS, tVt

所以理论上，如果我创建一些随机矩阵：

tA = np.array([
    [1, 0, 1],
    [0, 1, 0],
    [0, 1, 1],
    [0, 1, 0],
    [1, 1, 0]
])

然后调用

SVD

：

tU, tS, tVt = SVD(tA)

我应该能够将这些矩阵相乘并得到 tA，但是当我这样做时......

> print(tU @ tS @ tVt)
[[-3.33333333e-01  1.33333333e+00 -3.33333333e-01]
 [ 6.66666667e-01  3.33333333e-01  6.66666667e-01]
 [ 1.00000000e+00  1.00000000e+00  3.56682708e-16]
 [ 6.66666667e-01  3.33333333e-01  6.66666667e-01]
 [ 2.64708821e-16  1.00000000e+00  1.00000000e+00]]

只是为了验证这个 does 与 numpy 的 svd 方法一起工作：

tU, tSingularValues, tVt = np.linalg.svd(tA)
tS = np.zeros(tA.shape)
np.fill_diagonal(tS, tSingularValues)
print(tU @ tS @ tVt)

输出：

[[ 1.00000000e+00 -1.11022302e-16  1.00000000e+00]
 [-5.55652131e-19  1.00000000e+00  1.11577955e-16]
 [ 1.48107367e-16  1.00000000e+00  1.00000000e+00]
 [ 8.93209812e-17  1.00000000e+00  2.43745926e-16]
 [ 1.00000000e+00  1.00000000e+00  6.17126315e-17]]

你可以看到确实四舍五入到正确的输出。

我很确定我的 [S] 和 [V]^T 是正确的，但是我的方法中的 [U] 是关闭的；我的方法是为 [U]:

输出这个

[[-3.65148372e-01 -8.16496581e-01 -7.08401324e-17  6.05578786e-02
   4.47213595e-01]
 [-3.65148372e-01  4.08248290e-01  2.16496693e-16 -6.40036054e-01
   4.47213595e-01]
 [-5.47722558e-01 -1.93343331e-17 -7.07106781e-01 -6.05578786e-02
  -4.47213595e-01]
 [-3.65148372e-01  4.08248290e-01 -1.06808792e-16  7.61151811e-01
   4.47213595e-01]
 [-5.47722558e-01  6.39305215e-17  7.07106781e-01 -6.05578786e-02
  -4.47213595e-01]]

将其与 numpy 的 svd 进行比较[U]:

[[-3.65148372e-01  8.16496581e-01  5.67748493e-16  1.18391207e-01
  -4.31258069e-01]
 [-3.65148372e-01 -4.08248290e-01 -3.91737304e-16 -5.63487672e-01
  -6.18451004e-01]
 [-5.47722558e-01 -2.49196703e-16  7.07106781e-01 -1.18391207e-01
   4.31258069e-01]
 [-3.65148372e-01 -4.08248290e-01 -3.61832812e-16  8.00270086e-01
  -2.44065135e-01]
 [-5.47722558e-01  5.90601783e-16 -7.07106781e-01 -1.18391207e-01
   4.31258069e-01]]

你可以看到它们非常相似 - 除了其中一些莫名其妙地翻转了他们的标志。本质上，列向量 #2 和 #3 是错误的符号（无论如何舍入为 0 的值除外）。

我可以手动翻转这些列的标志：

tU = tU.T
tU[1] = -tU[1]
tU[2] = -tU[2]
tU = tU.T

现在输出基本上是正确的：

> print(tS @ tU @ tVt)
[[ 1.00000000e+00 -6.66133815e-16  1.00000000e+00]
 [ 4.30642036e-16  1.00000000e+00 -1.53086279e-16]
 [ 5.00222690e-16  1.00000000e+00  1.00000000e+00]
 [ 3.13052837e-16  1.00000000e+00  1.86547524e-16]
 [ 1.00000000e+00  1.00000000e+00  1.21168839e-16]]

但是如果不依赖 numpy 的 svd，我不知道我怎么能够提前知道哪些列向量需要反映出来，因为 numpy 将特征向量交给了我错误的符号。

有没有办法强制

np.linalg.eig

返回正确签名的向量？