特異値分解と固有値分解が同じ結果を出していることの確認（Rのsvd、eigen関数）

ある行列A（サンプルN×変数 D）があったとして、
$A = U \Sigma V^T$
のように行列分解するのが、特異値分解（SVD）。

Aの共分散行列
$S = \frac{1}{N-1}A^T A$
に対して、
$SV = V\Lambda$
となるような、固有ベクトル$V$を探す、すなわち、
$S = V\Lambda V^T$
のように行列分解するのが、固有値分解（EVD、Vは特異値分解のものと同じ）。

$S = \frac{1}{N-1}A^T A = \frac{1}{N-1}(U \Sigma V^T)^T U \Sigma V^T = \frac{1}{N-1}V \Sigma U^T U \Sigma V^T = \frac{1}{N-1}V \Sigma^2 V^T$
となることから、

SVDの特異値とEVDの固有値には以下のような関係があるはず。
$\frac{1}{N-1}\Sigma^2 = \Lambda$

というのは、数式上ではわかるが、実際にRの関数が同じ値を返すか不安だったので、確認してみた。

準備

まずデータを用意する。

# 10サンプル × 15次元のデータとする
set.seed(123)
A <- matrix(runif(150), nrow=10)
A <- scale(A, scale=FALSE) # あらかじめmean-centeringはしておく

続いて、SVD、EVDを実行

# SVD
res.svd <- svd(A)
# EVD
S <- t(A) %*% A / (10 - 1)
res.evd <- eigen(S)

固有値は同じか

まずは、SVDが返す特異値の自乗$\Sigma^2$と、EVDが返す固有値$\Lambda$が同じかを調べた

RではSVDはsvd関数、EVDはeigen関数があるので、これらを使った

# SVDの特異値^2 = EVDの固有値
plot(1:10, res.svd$d^2/(10 - 1), col="black", type="l", ylim=c(0,4), xlim=c(1,10), main="Eigen Value")
par(new=TRUE)
plot(1:10, res.evd$values[1:10], col="blue", type="l", ylim=c(0,4), xlim=c(1,10), ann=FALSE)
legend("topright", legend=c("SVD", "EVD"), col=c("black", "blue"), pch=16)

かなり一致しているのが見える

固有ベクトルは同じか

次に、固有ベクトル$V$がsvdとeigen関数で同じかを調べた

固有ベクトルは長さ1のベクトルなので、同じ方向を向いているほど、内積が1に近づく

180度反対の向きのものも正解なので、内積の絶対値を確認した

evec.svd <- res.svd$v
evec.evd <- res.evd$vectors[,1:10]
image(abs(t(evec.svd) %*% evec.evd), xlab="SVD", ylab="EVD", main="Eigen Vector")

最後が怪しいが、それ以外は一致した

ローディングは同じか

最後に、長さの違うもう一つの固有ベクトル$U$が同じになるかを調べた

SVDでは、
$A = U \Sigma V^T$
の$U$ですでに計算してある。

EVDでは、データ行列$A$と固有ベクトル$V$と、固有値$\Lambda$しかないが、$A$に$V$をかければ、$U$がもとまる

$A V = U \Sigma$

ただし、$\Sigma$分の重さがかかってるので、長さを1にして$\Sigma$を取り除く

loading.evd <- A %*% evec.evd
loading.svd <- res.svd$u
apply(loading.evd, 2, function(x){norm(x, "2")}) # この固有ベクトルにはまだ、固有値分の重みがかかっている
apply(loading.svd, 2, function(x){norm(x, "2")}) # SVDのVはΣと分離されているからOk
for(i in 1:ncol(loading.evd)){
	loading.evd[, i] <- loading.evd[, i] / norm(loading.evd[, i], "2")
}
image(abs(t(loading.evd) %*% loading.svd), xlab="EVD", ylab="SVD", main="Loading")

やはり、最後（10番目の固有ベクトル）が何かおかしい

固有ベクトルなら、互いに直行しているため、内積が0になる

そのため、SVDの$U$と、EVDの$U$で各々内積をみた

image(crossprod(loading.svd), main="SVD")
image(crossprod(loading.evd), main="EVD")

これを見る限り、EVDの10番目の軸が他の軸と直行していない

このように、原因はよくわからないが、eigen関数の高次の部分が不安定になる事例が見受けられた