いつもながらですが,R でプログラムを書くのは,非常に簡単。びっくりするほど。
gendat <- function(nc, r)
{
nv <- ncol(r)
# 仮のデータ行列を作る。この時点では変数間の相関は近似的に0である。
z <- matrix(rnorm(nv*nc), ncol=nv)
# 主成分分析を行い,主成分得点を求める。この時点で変数間の相関は完全に0である。
res <- eigen(r2 <- cor(z))
coeff <- solve(r2) %*% (sqrt(matrix(res$values, nv, nv, byrow=TRUE))*res$vectors)
z <- t((t(z)-apply(z, 2, mean))/sqrt(apply(z, 2, var)*(nc-1)/nc)) %*% coeff
# コレスキー分解の結果をもとに,指定された相関係数行列 を持つように主成分得点を変換する。
z %*% chol(r)
}
nc <- 20
z <- gendat(nc, matrix(c(1, 0.512, 0.335, 0.512, 1, 0.467, 0.335, 0.467, 1), ncol=3))
z
apply(z, 2, mean)
sqrt(apply(z, 2, var)*(nc-1)/nc)
cor(z)