多次元アレイ

  • こちらで多次元データを扱った。
  • そこでは、k次元立方格子状のデータ格納をarray()を使って行うとともに
  • それをハンドリングするためのn進数ハンドリングを考えた。
  • 多項分布では項数-1の次元で格納することができるはず・・・
  • お試し中ソース
k<-3 # カテゴリ数
# 自由変数の数はk-1
n<-10 # 総数
library(MCMCpack)

prob<-rdirichlet(1,rep(1,k)) # カテゴリの生起確率

aaa<-array(0,rep((n+1),(k-1))) # 計算結果格納用アレイ
# (n+1)^自由変数の数 の値がとりうる値
keta<-0:(k-2) # n進数計算のための道具
keta<-(n+1)^keta # n進数のアレイ内番地計算用の道具

counter<-rep(0,k-1) # 同上
counter[1]<--1
loop<-TRUE
while(loop){
	counter[1]<-counter[1]+1 #順々にn進数をその一桁目から増やす
	# くらい上がりの処理
	for(i in 1:(k-2)){
		if(counter[i]==n+2){
			counter[i]<-0
			counter[i+1]<-counter[i+1]+1
		}
	}
	# 1列のベクトルで言うと、これは何番目?
	banchi<-sum(keta*counter)+1
	vvv<-counter+1
	# vvvはそれぞれの自由変数の値
	# 不自由な変数の値は
	lastval<-n-sum(vvv)
	vvv2<-sum(vvv)-k+1
	tmpval<-c(vvv,vvv2)
	tmpprob<-prob^tmpval
	aaa[banchi]<-prod(tmpprob)*exp(lgamma(n+1)-sum(lgamma(tmpval+1)))
	
	if(sum(counter)==(n+1)*(k-1)){
		loop<-FALSE
	}
}