樹状分布

• こちらで濃淡がいろいろな分布とその上のディリクレ分布のことをやっている
• こちらでフローサイトメータをやっている
• まずまずな樹状形の分布をランダム発生させてみる
• 適当に枝分かれルールを作り、枝の長さをつくり、それを順次「生やす」ことで木型を作る
• 枝は生えるところから減衰するように粗密を作る
• 作成した木は梢に行くにしたがって薄くなる折れ線木
• これに乱雑項を入れる

library(vegan)
library(GPArotation)
library(igraph)
library(rgl)
library(MCMCpack)
# d:次元
# n.steps:枝分かれ回数
# n.pt.raneg:小集団の点の数
# NN:数を増やすスケール
# min.v,max.v:枝の長さの上下限
# できた木の上の点に誤差項を入れるときの係数
my.FACS.pts <- function(d,n.steps=10,n.pt.range=100:1000,NN=100,min.v=5,max.v=10,v=0.3){
	# 点の座標
	X <- matrix(0,0,d)
	# 枝分かれを繰り返す
	# ステップごとの点の数
	tmp.ns <- sample(n.pt.range,n.steps,replace=TRUE)
	# ステップごとの小集団の広がり(半径)を決める係数
	tmp.Ls <- runif(n.steps,min.v,max.v)
	for(i in 1:n.steps){
		tmp.n <- tmp.ns[i]
		tmp.X <- matrix(0,tmp.n,d)
		tmp.L <- tmp.Ls[i]
		# 基準点からの距離の２乗に応じて確率が下がるような生起確率を作り
		tmp.V <- seq(from=0,to=1,length=tmp.n*NN)
		tmp.prob <- (tmp.V-mean(tmp.V))^2
		tmp.prob <- tmp.prob/sum(tmp.prob)
		# それに応じてtmp.n個の点を発生させる
		tmp.S <- tmp.V[sort(sample(1:length(tmp.V),tmp.n,replace=TRUE,prob=tmp.prob))]*tmp.L
		# 枝上の位置を算出したので、それを適当な方向に回転する
		tmp.X[,1] <- tmp.S
		tmp.R <- Random.Start(d)
		tmp.X <- tmp.X %*% tmp.R
		# 枝は既存の位置から生やす
		if(i > 1){
			# 既存の位置から１点を選ぶ
			tmp.st <- sample(1:length(X[,1]),1)
			#tmp.st.2 <- sample(1:tmp.n,1)
			# 新しい枝は始点から「生える」だけにする
			tmp.st.2 <- 1
			st.pt <- X[tmp.st,]
			connect.pt <- tmp.X[tmp.st.2,]
			# 生え際を基準に座標を取り直す
			tmp.X <- t(t(tmp.X) + st.pt - connect.pt)
		}
		# 新枝の点を加える
		X <- rbind(X,tmp.X)
	}
	# 
	X <- X + rnorm(length(X),0,v)
	X.curve <- X
	for(i in 1:d){
		X.curve[,i] <- sign(X.curve[,i]) * abs(X.curve[,i]^2)
	}
	X <- X.curve
	return(X)
}
# 適当に線状の点の集まりを作る
d <- 3

n.steps <- 10
n.pt.range <- 100:1000
NN <- 100

my.X <- my.FACS.pts(d,n.steps,n.pt.range,NN)
plot(as.data.frame(my.X))