用numpy的eigh和svd计算的特征向量不匹配

考虑奇异值分解M = USV *。然后M * M的特征值分解给出M * M = V（S * S）V * = VS * U * USV *。我希望通过展示由eigh函数返回的特征向量是相同的svd功能恢复验证与numpy的这种平等：用numpy的eigh和svd计算的特征向量不匹配

import numpy as np 
np.random.seed(42) 
# create mean centered data 
A=np.random.randn(50,20) 
M= A-np.array(A.mean(0),ndmin=2) 

# svd 
U1,S1,V1=np.linalg.svd(M) 
S1=np.square(S1) 
V1=V1.T 

# eig 
S2,V2=np.linalg.eigh(np.dot(M.T,M)) 
indx=np.argsort(S2)[::-1] 
S2=S2[indx] 
V2=V2[:,indx] 

# both Vs are in orthonormal form 
assert np.all(np.isclose(np.linalg.norm(V1,axis=1), np.ones(V1.shape[0]))) 
assert np.all(np.isclose(np.linalg.norm(V1,axis=0), np.ones(V1.shape[1]))) 
assert np.all(np.isclose(np.linalg.norm(V2,axis=1), np.ones(V2.shape[0]))) 
assert np.all(np.isclose(np.linalg.norm(V2,axis=0), np.ones(V2.shape[1]))) 

assert np.all(np.isclose(S1,S2)) 
assert np.all(np.isclose(V1,V2))

最后断言失败。为什么？

来源

2015-01-05 matus

您可以添加对于所有对角线元素都是正数，即令M2 = M + a * I，其中a足够大以使M2正半边。那么SVD和eigh应该更好地达成一致。 –

只需使用小数字来调试您的问题。

开始A=np.random.randn(3,2)你有大小(50,20)

在我随机的情况下大得多的矩阵，而不是，我发现

v1 = array([[-0.33872745, 0.94088454], 
    [-0.94088454, -0.33872745]])

和v2：

v2 = array([[ 0.33872745, -0.94088454], 
    [ 0.94088454, 0.33872745]])

，他们只对不同一个标志，显然，即使标准化为单位模块，向量也可能因符号而不同。

现在，如果你尝试的伎俩

assert np.all(np.isclose(V1,-1*V2))

为原来的大矩阵，失败......一次，这是确定。会发生什么是一些向量乘以-1，其他一些则没有。

检查向量之间平等的一个正确的方法是：

assert allclose(abs((V1*V2).sum(0)),1.)

而事实上，得到的是如何工作的，你的感觉可以打印此量：

(V1*V2).sum(0)

那的确是根据载体：+1或-1：

array([ 1., -1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 
    1., -1., 1., 1., 1., -1., -1.])

编辑：这将发生在大多数情况下，特别是如果从随机矩阵开始。但是请注意，本次测试将有可能，如果失败的一个或多个特征值具有尺寸比1大的本征空间，如下面他的评论中指出@Sven Marnach：

有可能不仅仅是矢量乘以其他方面的差异-1。如果任何特征值具有多维特征空间，你可能会得到固有空间的任意正交基，并这样的基地可能对彼此的arbitraty 一神教矩阵旋转

来源

2015-01-05 15:08:13 gg349

@matus好的，我迷路了:)但我相信你的判断，所以我会删除我的意见，不要混淆未来的读者。干杯! – BartoszKP

除了向量乘以-1之外，可能还有其他差异。如果任何特征值具有多维特征空间，那么可以得到该特征空间的任意正交基，并且这样的基可以通过任意一元矩阵相互旋转。 –

@SvenMarnach，这是一个非常有用的观点。我将编辑这篇文章来指出这个警告 – gg349

用numpy的eigh和svd计算的特征向量不匹配

回答

相关问题