2011-08-29

条件付きself-avoiding path

Self-avoiding pathは「すでに通った地点」は通過してはいけないというルールでのランダムウォークのpath(こちら)
モジホコリ(こちら)は、自身のすでに通った地点を避ける(自身の排せつ物があるから？餌があるとは考えにくいから？)
感染症(こちら)では、自身が通過した地点(ホスト)は免疫力を有するがために、再度は通過できない(永久の免疫力なら)
生物がどんなに昔の過去にまで遡って情報を使用することはない(はず)で、免疫力も低下する(特に、繰り返し感染するような感染症の場合には)
Self-avoiding path的には、「ある一定の過去」までの記録上、通過したことになっているときににそこを避ける
双曲幾何的なサンゴ(こちら)の場合は、「通過した点には『自分』が居る」からそこには入れない、という制約があって、その上で、サンゴ虫が増殖するので、「存在場所を求めて」広がっていく
膜活動電位を用いる系も不応期という「self-avoiding」な条件を持つ(こちら(心筋電気伝達)やこちら(神経生理))
どれも「条件付きself-avoiding walk」
ひとまず、一番単純な「２次元」「１匹」「増殖なし」「一番昔まで記憶を辿る」という条件でのSelf-avoiding生物をプロットしてみる
- 行き場のないどん詰まりも生じる条件

L<-50
k<-2
B<-array(0,rep(2*L+1,k))

zero<-L+1

B[zero,zero]<-2

r<-0.5
p<-c(0,r,1)

T<-100
par(ask=FALSE)
for(i in 2:T){
	image(B,col=gray(c(2,1.5,0)/2))
	Sys.sleep(0.05)
	newB<-B
	exist<-which(B==2,arr.ind=TRUE)
	#print(exist)
	s<-sample(1:length(exist[,1]))
	#print(s)
	for(j in 1:length(exist[,1])){
		tmpr<-runif(1)
		n<-0
		for(j2 in 1:length(p)){
			if(tmpr<p[j2]){
				n<-j2
			}
		}
		ok<-matrix(rep(exist[s[j],],2*k),ncol=2*k)
		#print(ok)
		okid<-c()
		for(j2 in 1:k){
			ok[j2,2*(j2-1)+1]<-ok[j2,2*(j2-1)+1]-1
			ok[j2,2*j2]<-ok[j2,2*j2]+1
			if(newB[1+sum((ok[,2*(j2-1)+1]-1)*(2*L+1)^(0:(k-1)))]==0){
				okid<-c(okid,2*(j2-1)+1)
			}
			if(newB[1+sum((ok[,2*j2]-1)*(2*L+1)^(0:(k-1)))]==0){
				okid<-c(okid,2*j2)
			}
		}
		if(length(okid)>0){
			here<-sample(okid,1)
			newB[1+sum((ok[,here]-1)*(2*L+1)^(0:(k-1)))]<-2
			newB[1+sum((exist[s[j],]-1)*(2*L+1)^(0:(k-1)))]<-1
		}
		
	}
	B<-newB
	
}

どん詰まりを回避するとすると次のような条件が簡単に思いつく
- 後退してよいという条件
- 遡る過去を限定するという条件
過去を限定してみよう
様子をanimationパッケージのsaveGIF()関数(こちら)で映画にしよう

# 過去をある程度忘れる
library(animation)
L<-30
k<-2
B<-array(0,rep(2*L+1,k))

zero<-L+1

active<-50
B[zero,zero]<-active

r<-0.5
p<-c(0,r,1)

T<-500
par(ask=FALSE)
saveGIF({
for(i in 2:T){
	print(i)
	Sys.sleep(0.01)
	newB<-B-1
	newB[which(newB<0)]<-0
	newB[which(newB==active-1)]<-active
	image(newB,col=gray((active:0)/active))
	exist<-which(B==active,arr.ind=TRUE)
	print(exist)
	s<-sample(1:length(exist[,1]))
	#print(s)
	for(j in 1:length(exist[,1])){
		tmpr<-runif(1)
		n<-0
		for(j2 in 1:length(p)){
			if(tmpr<p[j2]){
				n<-j2
			}
		}
		ok<-matrix(rep(exist[s[j],],2*k),ncol=2*k)
		print(ok)
		okid<-c()
		for(j2 in 1:k){
			ok[j2,2*(j2-1)+1]<-ok[j2,2*(j2-1)+1]-1
			ok[j2,2*j2]<-ok[j2,2*j2]+1
			print(ok[,2*(j2-1)+1])
			print(ok[,2*j2])
			if(newB[1+sum((ok[,2*(j2-1)+1]-1)*(2*L+1)^(0:(k-1)))]==0){
				okid<-c(okid,2*(j2-1)+1)
			}
			if(newB[1+sum((ok[,2*j2]-1)*(2*L+1)^(0:(k-1)))]==0){
				okid<-c(okid,2*j2)
			}
		}
		if(length(okid)>0){
			here<-sample(okid,1)
			newB[1+sum((ok[,here]-1)*(2*L+1)^(0:(k-1)))]<-active
			newB[1+sum((exist[s[j],]-1)*(2*L+1)^(0:(k-1)))]<-active-1
		}
		
	}
	B<-newB
	
}
},interval=0.05)

少し簡単にして、cpp化してRcppでつないでみる→こちら

# k:次元,past:考慮する過去(-1は無限に),t:歩数,max.t:最大トライ数
selfAvoidingWalkR <- function(k,past=-1,t=100,max.t=1000){
	X <- matrix(0,max.t,k)
	R <- sample(1:k,max.t,replace=TRUE)
	R.sign <- sample(c(-1,1),max.t,replace=TRUE)
	cnt <- 1
	for(i in 1:(max.t-1)){
		here <- X[cnt,]
		toward <- X[cnt,]
		toward[R[i]] <- toward[R[i]] + R.sign[i]
		if(past == 0){
			cnt <- cnt+1
			X <- X[cnt,] <- toward
		}else{
			if(past == -1){
				tmp.past <- 1
			}else{
				tmp.past <- cnt-past+1
				tmp.past <- max(1,tmp.past)
			}
			if(cnt==tmp.past){
				tmp.X <- matrix(X[cnt,],nrow=1)
			}else{
				tmp.X <- X[tmp.past:cnt,]
			}
			tmp.d <- apply(abs(t(tmp.X)-toward),2,sum)
			if(prod(tmp.d)!=0){
				cnt <- cnt+1
				X[cnt,] <- toward
			}
		}
		if(cnt==t){
			break
		}
	}
	X[1:cnt,]
}


k<-3
max.t<-10000
past <- -1
t <- 1000
B <- selfAvoidingWalkR(k,past=past,t=t,max.t=max.t)
plot(B,type="b")