面白いなぁ。

barnacle<-function(b0,a0,a1,b1,additive_var){
z<-c()
z[1]=0.1
maxt=1000
for(t in 1:(maxt-1)){
W<-z[t]*(z[t]/2*b1+(1-z[t]/2)*b0)+(1-z[t])**1
z[t+1]<-z[t]+deltaz
}
plot(z,type="l",xlab="generation",ylab="swichpoint",ylim=c(0,1))
text(500,0.1,paste("z=",round(z[1000],4),"equ=",round*2,cex=1.1)
}

*1:1-z[t])/2*a1+(1-z[t])/2*a0) deltaz<-additive_var/W*(z[t]*(b1-a1)+(1-z[t])*(b0-a0

*2:a0-b0)/(b1-b0+a0-a1),4