2013-02-23

決断関数の構造、カテゴリカルな決断

上の例では、複数の選択肢から発生するすべての帰結に順序を入れた。さらに、その順序に応じて、なるべく、「順序が前」の帰結が得られるように「単純な」処理を繰り返した
「順序〜大小関係」のみを用いるというのは、結局のところ、不等号には重み１の差をつけ、等号には重みの差をつけずに、「線形な面」をつくっているのと(多分)同じ
そうすると、「治療選択における、治癒と死亡」ただし、「治癒にも副作用ありとなしがあって、死亡にも副作用ありなしがある」のような場合には、「副作用の有無によらず治癒」は「副作用の有無によらず死亡」の上位におきたい
すべての帰結を同一尺度に並べ、ただし、不等号に持たせる重みに軽重をつける、というのも一つのやり方だろう
また別のやり方としては、「順序〜大小関係」にヒエラルキーを入れることになるだろう
いずれにしろ、その「関数・重みづけ係数」がわかっていれば、簡単なわけで、それがわからないながら、前進したい、と、それをする方便を入れないといけない
そんな方向としては、モンテカルロ比較をしながら、「選択肢Aだと、こうなって、選択肢Bだとこうなる、と仮定したら、AとBとどっちをとりますか？」だけでなく(これだと順序を入れただけ)、どのくらい強く(主観的で酔い)Aをとりたいですか？」とインターラクティブにモンテカルロを回す、というのもありかもしれない。
ひとまず、２治療、副作用持ち込み編(順序だけで評価)の書き散らしソースを

# オプション数
k <- 2
# 各オプションの帰結の種類数
s <- sample(2:5,k,replace=TRUE)
# 各オプション３個にしてみる
s <- rep(4,k)
s.list <- list()
for(i in 1:k){
	s.list[[i]] <- 1:s[i]
}
# 各オプションごとに帰結の生起確率
true.prob <- list()
for(i in 1:k){
	true.prob[[i]] <- rdirichlet(1,rep(10,s[i]))
}
# ちょっと固定してみる
# オプション２とオプション３とが競り合うように
#true.prob[[1]] <- c(0.5,0.5)
#true.prob[[2]] <- c(0.45,0.55)
#true.prob[[3]] <- c(0.5,0.1,0.4)
#true.prob[[2]] <- c(0.0,0.5,0.5)
true.prob[[1]] <- c(0.0,0.6,0.1,0.3)
true.prob[[2]] <- c(0.4,0.15,0.15,0.3)
# 起こりうる帰結の組合せの網羅
cat.combin <- expand.grid(s.list)

# どの帰結を好むかの行列を作る
#preference <- sample(1:k,prod(s),replace=TRUE)
preference <- matrix(1,prod(s),k)
for(i in 1:prod(s)){
	tmp.s <- sample(0:(k-1),1)
	preference[i,sample(1:k,tmp.s)] <- 0
	preference[i,] <- preference[i,]/sum(preference[i,])
}
#preference <- matrix(c(1,1,0,1,1,0,1,1),ncol=2)
# ホタルのモデル
scores <- matrix(c(7,4,3,1,6,5,3,2),ncol=2)
for(i in 1:length(preference[,1])){
	tmp <- c(scores[cat.combin[i,1],1],scores[cat.combin[i,2],2])
	preference[i,] <- as.numeric(tmp==max(tmp))
}
# 和が１となるように標準化
for(i in 1:length(preference[,1])){
	preference[i,] <- preference[i,]/sum(preference[i,])
}
# 生起確率は全く不明なところからスタート
pre <- list()
for(i in 1:k){
	pre[[i]] <- sample(0:0,s[i],replace=TRUE)
}
# そのつど、確率的によさげなオプションを選んでは、その結果をhistoryに記録し
# 記録に応じて、次の選択を判断する
n.history <- 1000
history <- matrix(0,n.history,sum(s))
# 結果からディリクレ乱数で「生起確率」を推定する
# 作成するディリクレ乱数の個数
n.iter <- 100
#N.iter <- 1

for(h in 1:n.history){
	tmp <- c()
	for(j in 1:k){
			tmp <- c(tmp,pre[[j]])
	}
	history[h,] <- tmp
	chosen <- rep(0,k)
	#for(ii in 1:N.iter){
		this.time.chosen <- rep(0,k)
		rs <- list()
		for(j in 1:k){
			rs[[j]] <- rdirichlet(n.iter,pre[[j]]+1)
		}

		# 推定した生起確率とディリクレ乱数から
		# もっとも好みに沿った結果が得られるオプションを
		# 確率的に選ぶ
		for(i in 1:n.iter){
			tmp.rs <- list()
			for(j in 1:k){
				tmp.rs[[j]] <- rs[[j]][i,]
			}
			probs <- expand.grid(tmp.rs)
			tmp.prob <- apply(probs,1,prod)
			tmp2 <- preference * tmp.prob
			tmp3 <- apply(tmp2,2,sum)
			tmp.max <- which(tmp3 == max(tmp3))
			if(length(tmp.max)==1){
				choice <- tmp.max
			}else{
				choice <- sample(tmp.max,1)
			}
			#choice <- tmp.max[sample(seq(tmp.max),1)]
			this.time.chosen[choice] <- this.time.chosen[choice]+1
			#print(choice)
			#tmp <- matrix(NA,N.iter,k)
			#for(j in 1:k){
			#	tmp[,j] <- sample(1:s[j],replace=TRUE,prob=rs[[j]])
			#}
		}
		tmp.max <- which(this.time.chosen==max(this.time.chosen))
		if(length(tmp.max)==1){
			selected.chosen <- tmp.max
		}else{
			selected.chosen <- sample(tmp.max,1)
		}
		selected.chosen <- tmp.max[sample(seq(tmp.max),1)]
		chosen[selected.chosen] <- chosen[selected.chosen]+1
	#}
	#tmp.max <- which(chosen==max(chosen))
	#if(length(tmp.max)==1){
	#	really.chosen <- tmp.max
	#}else{
	#	really.chosen <- sample(tmp.max,1)
	#}
	# 選んだオプションの結果は、隠された真の生起確率によって観察される
	really.chosen <- sample(1:k,1,prob=this.time.chosen/sum(this.time.chosen))
	#really.chosen <- tmp.max[sample(seq(tmp.max),1)]
	really.happen <- sample(1:s[really.chosen],1,prob=true.prob[[really.chosen]])
	pre[[really.chosen]][really.happen] <- pre[[really.chosen]][really.happen]+1

}

matplot(history,type="l")