クオンタイル

  • 乱数を使う
my.qdirichlet <- function(v,a,n=10000){
	r <- rdirichlet(n,a)
	d <- ddirichlet(r,a)
	d0 <- ddirichlet(v,a)
	length(which(d<=d0))/n
}
  • ベータ分布で試してみる

n <- 100
a <- c(1,1)
a2 <- c(3,8)
R <- rdirichlet(n,a)
q <- rep(0,n)
for(i in 1:length(R[,1])){
	q[i] <- my.qdirichlet(R[i,],a2,n=100)
}
plot(pbeta(R[,1],a2[1],a2[2],lower.tail=FALSE),q)
  • 3項で試してみる

n <- 1000
a <- c(1,1,1)
a2 <- c(3,8,10)
R <- rdirichlet(n,a)
q <- rep(0,n)
for(i in 1:length(R[,1])){
	q[i] <- my.qdirichlet(R[i,],a2,n=100)
}

M <- matrix(c(1,0,-1/2,sqrt(3)/2,-1/2,-sqrt(3)/2),2,3)
R. <- R %*% t(M)
Ctr <- (a2/sum(a2)) %*% t(M)
plot(R.,col=gray(1-q),pch=20)
points(Ctr,pch=20,col=2)