特異値分解・固有値分解系の分析

主成分分析とか因子分析とかの手法は特異値分解系でまとめられますが、調べてみると意外に良い資料がなくてなかなか理解が進まないですorz しかしとりあえず調べれる範囲で理解した部分を記事にしておきます。

教科書はこちら。私が学部4年のときに読んだ本で、自分にとっては一番分かり易く書かれています。しかし絶版orz

多次元データの解析 (シリーズ入門統計的方法 3)

多次元データの解析 (シリーズ入門統計的方法 3)


※以下記事ですが、まだ曖昧な部分も多いです。間違いがあればご指摘して下されば幸いです。「コレスポンデンス分析」と「多次元尺度法」はもう少し理解を深めるために、別途調査する予定です。



【特異値分解・固有値分解って？】

  • 特異値分解
    • 行列Aを次のように分解する。
    • A=UDV
    • Dは対角成分以外がゼロの行列で、対角成分を特異値と呼ぶ。UやVのベクトル成分を特異ベクトルと呼ぶ。
# データの共分散行列を特異値分解
svd(var(USArrests))

# データの共分散行列を固有値分解
eigen(var(USArrests))

# データを標準化せずに主成分分析
prcomp(USArrests)

どれも同じ結果になる（固有値分解は正方行列しか扱えない）。主成分分析はデータの共分散行列を特異値分解したものということが分かる。また次に示すように、標準化したデータの主成分分析は、相関行列の特異値分解である（標準化したデータの特異値分解とも一致する）。

svd(cor(USArrests))
svd(scale(USArrests))
prcomp(USArrests, scale = TRUE)

主成分分析は固有値分解ないしは特異値分解である。




【記述的な分析】

prcomp(USArrests, scale = TRUE)
summary(prcomp(USArrests, scale = TRUE))
biplot(prcomp(USArrests, scale = TRUE))

  • 因子分析
    • 潜在変数を仮定して次元縮小
      • x1=a11f1+a12f2+…+e1
      • x2=a21f1+a22f2+…+e2
    • f：共通因子、e：特殊因子（誤差）
    • 標準化した行列の共分散行列Sを次のように分解する
      • S=AA’+ε
    • Aは行列Xに共通因子fがあると仮定した場合の係数行列、εは誤差（特殊因子）の二乗を対角成分とする対角行列である
    • 数学的な問題点：分散であるeが負になる（不適解）場合がある→条件を付けて最適化。Aが一意でない→これを利用して、解釈しやすいように「回転」させる。反復法によって推定するため、アルゴリズムや初期値によってはパラメータが変化する可能性がある。
    • 回転させ方：1つの変数は1つの共通因子に対してだけ、負荷量（相関）が大きい。それぞれの共通因子の負荷量は、-1か1または0に近く、中間値をなるべくとらない。
    • 共通因子同士が独立（無相関）になるように回転させる→直交回転（ヴァリマックス法）。相関があっても良い→斜交回転（プロマックス法）。
v1 <- c(1,1,1,1,1,1,1,1,1,1,3,3,3,3,3,4,5,6)
v2 <- c(1,2,1,1,1,1,2,1,2,1,3,4,3,3,3,4,6,5)
v3 <- c(3,3,3,3,3,1,1,1,1,1,1,1,1,1,1,5,4,6)
v4 <- c(3,3,4,3,3,1,1,2,1,1,1,1,2,1,1,5,6,4)
v5 <- c(1,1,1,1,1,3,3,3,3,3,1,1,1,1,1,6,4,5)
v6 <- c(1,1,1,2,1,3,3,3,4,3,1,1,1,2,1,6,5,4)
m1 <- cbind(v1,v2,v3,v4,v5,v6)
a <- factanal(m1, factors = 3) # varimax is the default
a$loadings
factanal(m1, factors = 3, rotation = "promax")

# 主成分分析と比較
# 因子分析と主成分分析は要約され方が違う→回転させてるから
prcomp(m1)
summary(prcomp(m1))

prcomp(m1, scale=T)
summary(prcomp(m1, scale=T))

  • コレスポンデンス分析
    • 分割表の正準相関分析として研究が始まる？（1930年代～、Fisher）
    • 縮約する次元数を指定する（2次元に縮約するのが一般的）
    • 数量化III類、双対尺度法（dual scaling）
    • 1960年代あたりからフランスでanalyse des correspondences（correspondence analysis）と命名→国際共通語に
    • 1974、M.O.Hillが「無視された多変量解析手法」として紹介
library(MASS)
a <- corresp(caith, nf=2)
biplot(a)


library(fastICA)
library(MASS)
x <- mvrnorm(n = 1000, mu = c(0, 0), Sigma = matrix(c(10, 3, 3, 1), 2, 2))
x1 <- mvrnorm(n = 1000, mu = c(-1, 2), Sigma = matrix(c(10, 3, 3, 1), 2, 2))
X <- rbind(x, x1)

a <- fastICA(X, 2, alg.typ = "deflation", fun = "logcosh", alpha = 1,
             method = "R", row.norm = FALSE, maxit = 200, 
             tol = 0.0001, verbose = TRUE)

par(mfrow = c(1, 3))
plot(a$X, main = "Pre-processed data")
plot(a$X%*%a$K, main = "PCA components")
plot(a$S, main = "ICA components")


loc <- cmdscale(eurodist)
x <- loc[, 1]
y <- -loc[, 2]

plot(x, y, type = "n", xlab = "", ylab = "", asp = 1, axes = FALSE,
     main = "cmdscale(eurodist)")
text(x, y, rownames(loc), cex = 0.6)




【予測的な分析】

  • PLS（partial least squares）
    • 主成分分析＋予測
    • 予測変数も多次元可能なはずだが、Rでは対応していない？（SASは対応していた）
    • 一般化PLSはgplsパッケージ
library(pls)
data(yarn)
yarn.pls <- plsr(density ~ NIR, 6, data = yarn, validation = "CV")

library(gpls)

  • 正準相関分析（canonical correlation analysis）
    • 行列Xと行列Yを縮約してそれぞれ相関の高い変数を作る
library(CCA)
data(nutrimouse)
X=as.matrix(nutrimouse$gene[,1:10])
Y=as.matrix(nutrimouse$lipid)
res.cc=cc(X,Y)
plot(res.cc$cor,type="b")
plt.cc(res.cc)

ページTOPへ