ICA, Independent Component Analysis（独立成分分析）は、多変量の信号を複数の加法的な成分に分離するための計算手法です。
大量の多変量データを解析するとき、主成分分析をよく使っています。多次元データをばらつきの大きさ順に直交座標に変換していく主成分分析の方法は直観的に理解しやすく、また、取り組むべき解析対象の優先付けをし易くしてくれるので、大変妥当な手法だと思っています。
それでも、もっと良い方法が無いか、あるいは改善できるところがないか検討を重ねているのですが、「独立成分分析」という語を目にした時に？が思い浮かびました。（線型）独立って直交していることでないのかと勘違いをしていたのです。たまたま読んでいたサイト [1] の筆者も最初はそう考えたようでした。
直交関係にある複数のベクトルは線型独立の関係にありますが、逆は必ずしも成り立たないのです。が、その説明は省略します。
[1] のサイトで紹介されていた「独立成分分析」の R パッケージ fastICA [2] を使ってみました。
三つのデモプログラム
まずは、マニュアルの 4 ページに紹介されているデモサンプルを試してみました。
最初の例は、二次元の一様乱数を線形変換したサンプルを主成分分析と独立成分分析を適用した例です。
Example 1: un-mixing two mixed independent uniforms
library(fastICA)
S <- matrix(runif(10000), 5000, 2)
A <- matrix(c(1, 1, -1, 3), 2, 2, byrow = TRUE)
X <- S %*% A
a <- fastICA(X, 2, alg.typ = "parallel", fun = "logcosh", alpha = 1,
method = "C", 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")
runif で一様乱数を 10,000 個生成して、それを 5,000 行 × 2 列の行列にして、2 × 2 の行列を掛けて線型変換（回転、拡大）したデータをサンプルとして使っています。
Example 1 の実行例
左が元のデータ、真ん中が主成分分析で変換した結果、そして右が独立成分分析で変換した結果です。真ん中の主成分分析では、直交成分である PC1 と PC2 に変換したからと言って、データがその軸に沿って並ぶとは限らないことが判ります。
右の独立成分分析では、元のデータ S の分布の成分と同等な成分に変換しています。
二番目のサンプルは、いわゆるブラインド信号源分離をする例です。
Example 2: un-mixing two independent signals
library(fastICA)
S <- cbind(sin((1:1000)/20), rep((((1:200)-100)/100), 5))
A <- matrix(c(0.291, 0.6557, -0.5439, 0.5572), 2, 2)
X <- S %*% A
a <- fastICA(X, 2, alg.typ = "parallel", fun = "logcosh", alpha = 1,
method = "R", row.norm = FALSE, maxit = 200,
tol = 0.0001, verbose = TRUE)
par(mfcol = c(2, 3))
plot(1:1000, S[,1 ], type = "l", main = "Original Signals", xlab = "", ylab = "")
plot(1:1000, S[,2 ], type = "l", xlab = "", ylab = "")
plot(1:1000, X[,1 ], type = "l", main = "Mixed Signals", xlab = "", ylab = "")
plot(1:1000, X[,2 ], type = "l", xlab = "", ylab = "")
plot(1:1000, a$S[,1 ], type = "l", main = "ICA source estimates", xlab = "", ylab = "")
plot(1:1000, a$S[, 2], type = "l", xlab = "", ylab = "")
オリジナルの信号は、正弦波と三角波のような波形を持った二次元の信号 S です（下左）。これに 2 x 2 行列 A を乗じて、それぞれの成分を混ぜて線形変換します。これを X として各成分をプロットしたのが下の真ん中です。
Example 2 の実行例
X に対して独立成分分析して信号を分離したのが上右のプロットです。符号は逆になっていますが、元の信号に分離されているのが確認できます。
三番目のサンプルは、混合した二変数の正規分布を独立成分分析して射影追跡した例です。射影追跡とは， 多次元データの特徴的な構造を、低次元空間に線形射影することによって探るデータ解析手法です [3] 。
Example 3: using FastICA to perform projection pursuit on a mixture of bivariate normal distributions
library(fastICA)
if(require(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")
}
二つの異る二次元の正規分布群 x と x1 を合わせたデータ（下左）に対して主成分分析（下の真ん中）と独立成分分析（下右）を行っています。
Example 3 の実行例
比較している内容は判るのですが、どうせなら二つの群を色分けして確認したかったので、少し手を加ました。
Example 3 (b):
library(fastICA)
library(MASS)
plot.xy <- function(m, title, labels) {
colnames(m) <- labels
df = data.frame(m)
df$Group <- sapply(strsplit(as.character(row.names(df)), "_"), "[[", 1)
plot(df[, 1:2], main = title, type="n")
colors <- data.frame(Group = unique(df$Group))
colors$Color = c("#FF4000", "#0040FF")
for(g in unique(df$Group)) {
points(df[df$Group == g, 1:2], col = colors[colors$Group == g, "Color"], cex = 0.25)
}
}
n.max <- 1000
x1 <- mvrnorm(n = n.max, mu = c(0, 0), Sigma = matrix(c(10, 3, 3, 1), 2, 2))
row.names(x1) <- paste("A", c(1:n.max), sep = "_")
x2 <- mvrnorm(n = n.max, mu = c(2, -2), Sigma = matrix(c(10, 3, 3, 1), 2, 2))
row.names(x2) <- paste("B", c(1:n.max), sep = "_")
X <- rbind(x1, x2)
colnames(X) = c("x1", "x2")
a <- fastICA(X, 2, alg.typ = "deflation", fun = "logcosh", alpha = 1,
method = "R", row.norm = FALSE, maxit = 200,
tol = 0.0001, verbose = TRUE)
b <- prcomp(X)
par(mfrow = c(1, 3))
plot.xy(X, "Raw data", c("x", "y"))
#plot.xy(a$X %*% a$K, "PCA components", c("PC1", "PC2"))
plot.xy(b$x, "PCA components", c("PC1", "PC2"))
plot.xy(a$S, "ICA components", c("IC1", "IC2"))
Example 3 (b) の実行例
二種類のデータを混合して、主成分分析と独立成分分析をして比較する限り、独立成分分析の方がデータ群の差異を説明するための解釈の助けになるように思いました。もっと次元が多い複雑な系でどうなるか評価して、取り組んでいるデータ解析に適用できるかどうか検討したいと思います。
参考 サイト
- RでPCAとICA
- CRAN - Package fastICA
- 射影追跡：その考え方と実際
にほんブログ村
