ということで，なんで双対尺度法のコードまで作っちゃったかというと，全情報分析がやりたかったから。
すでに掲載した私流の双対尺度法コードの返り値をつかって，全情報分析に使う距離行列を算出するのが次のコード。
TIA.dist <- function(dat){
ndim <- dat$ndims
nr <- ncol(dat$nr)
nc <- ncol(dat$nc)
rname <- colnames(dat$nr)
cname <- colnames(dat$nc)
y <- t(dat$nr)
x <- t(dat$nc)
eta2 <- dat$eig
eta <- dat$sing
dist.mat <- matrix(0,nrow=nr+nc,ncol=nr+nc)
colnames(dist.mat) <- c(rname,cname)
rownames(dist.mat) <- c(rname,cname)
for(d in 1:ndim){
for(i in 1:nr){
# x-x
for(j in i:nr){
dist.mat[j,i] <- dist.mat[j,i]+ (eta2[d]*((y[i,d]-y[j,d])^2))
}
# x-y
for(j in 1:nc){
dist.mat[nr+j,i] <- dist.mat[nr+j,i]+(eta2[d]*(y[i,d]^2+x[j,d]^2-(2*eta[d]*y[i,d]*x[j,d])))
}
}
# y-y
for(i in 1:nc){
for(j in i:nc){
dist.mat[nr+j,nr+i] <- dist.mat[nr+j,nr+i]+(eta2[d]*(x[i,d]-x[j,d])^2)
}
}
}
dist.all <- sqrt(dist.mat)
dist.row <- dist.all[1:nr,1:nr]
dist.col <- dist.all[(nr+1):(nr+nc),(nr+1):(nr+nc)]
dist.rc <- dist.all[(nr+1):(nr+nc),1:nr]
return(as.dist(dist.all,upper=T))
}
ちなみに全情報分析とは，双対尺度法で最適な重みをつけられた行・列の空間の距離を計算することで全ての次元における全ての最適重みをつかって行・列データの分類を行う方法で，この関数で算出された距離行列に対して例えばk-means法を使うことを提案している。
もちろん，k-means法に限らなくても，改良k-means法や自己組織化マップ，MDSなどを使っても良いだろう（因子分析など次元圧縮を行うのは，せっかく全部の情報を使おうとしているこの手法の哲学に沿わないだろうからふさわしくないと思われるが）。
テキストの例を次のように入れてみたら，きちんと結果は再現された。もっとも，k-means法は乱数アルゴリズムを使っているから，必ず同じ結果が再現されるとは限らないが。
Rorschaha <- matrix(c(33,10,0,7,2,5,0,3,2,0,1,
10,5,2,0,9,9,3,2,1,4,6,
18,2,1,13,30,1,4,6,4,2,1,
1,1,26,1,4,2,5,2,1,1,0,
2,0,5,4,1,1,5,2,18,2,1,
6,0,18,2,6,1,21,3,2,2,0),nrow=11)
rownames(Rorschaha) <- c("Bat","Blood","Butterfly","cave","cloud","fire","fell","mask","Mt","rock","fog")
colnames(Rorschaha) <- c("Fear","Anger","Depress","Ambit","secure","affair")
ret <- DualScaling(Rorschaha)
tia <-TIA.dist(ret)
kmeans(tia,3)
皆さんもぜひ遊んでみてください。
出典；行動科学のためのデータ解析―情報把握に適した方法の利用
追記 関数名に誤記がありましたので直しました。前記事のDualScaling関数の戻り値をTIA.dist関数に渡すと距離行列が返ってきます。