numpy.linalg.eigでSVD(特異値分解)を行う方法

Question

np.linalg.svdを使わずnp.linalg.eig　で特異値分解を行おうとしております。
eigとsvdで行った場合、特異ベクトルの符合が揃いません。
eigのみで汎用的に特異値分解ができるように実装するのは出来ないのでしょうか。

import numpy as np
A = np.array([[2,4,1,3],[1,5,3,2],[5,7,0,7]])

#比較のためSVDでの実装
u, s, v = np.linalg.svd(A,full_matrices=False) 

#svd:u - eig:u_vec
#svd:s - eig:sigma
#svd:v - eig:v_vec

#左特異ベクトルと固有値
vals,vec =  np.linalg.eig(A@A.T)
#エルミネート行列ではないがeighで行うとSVDと同じになる
idx = np.argsort(vals)[::-1]

sigma = np.diag(np.sqrt(vals[idx]))
sigma = np.nan_to_num(sigma,0)
s = np.diag(s)

print('EIG_∑\n',sigma)
print('SVD_∑\n',s)

#左特異ベクトル
u_vec = vec[:,idx]
print('\nEig_左特異ベクトル\n' , u_vec)
print('SVD_左特異ベクトル\n' , u)


#右特異ベクトル
vals,vec =  np.linalg.eig(A.T@A)
#エルミネート行列ではないがeighで行うとSVDと同じになる
idx = np.argsort(vals)[::-1]
v_vec = vec[:,idx].T
v_vec = v_vec[:3,:]
print('\nEig_右特異ベクトル\n' , v_vec)
print('SVD_右特異ベクトル\n' , v)

#元配列に戻す
print('\n元配列A:\n',A)
print('\nSVD:\n' + str(np.round(u@s@v,2)))
print('\nEIG:\n' + str(np.round(u_vec@sigma@v_vec,2)))

Answer 1

eighを使いましょう。

#左特異ベクトルと固有値
vals,vec =  np.linalg.eig(A@A.T)
#エルミネート行列ではないがeighで行うとSVDと同じになる

上記のコメントは勘違いされています。A@A.T （およびA.T@A）は、Aが実数の時、常に対角行列（エルミート行列）です。

eigとeighは、エルミート行列に対しては、速度やソートの違いだけかと思いきや、eigには精度の問題もあるようです。そのため、エルミート行列においては、eighを使うことが推奨されます。

参考: Numpy で行列を対角化する際の精度

Comments

[t05708hn]

ありがとうございます.eig・eighどちらにしても行列が大きくなると精度が出ないのですね。勉強になりました。