Saturday, March 5, 2011

BM simulation with bounds

Inspired by another R-sig-phylo listserve discussion, I just added the capacity to conduct bounded Brownian simulation to my fastBM() function. This can be found on my R-phylogenetics page (direct link to code here).

Doing this was easy. I just added a simple while() loop which runs after a descendant node phenotype is calculated. The while() condition evaluates whether the phenotype is on the interval given by bounds, and if not, then it returns the phenotype to within the bounds by bouncing it back by the amount that it would have otherwise exceeded the boundary.

while(y[i,2] < bounds[1] || y[i,2] > bounds[2]){
     if(y[i,2]< bounds[1]) y[i,2] <- 2*bounds[1]-y[i,2]
     if(y[i,2] > bounds[2]) y[i,2] <- 2*bounds[2]-y[i,2]

Seems to work . . . . Feedback welcome!

No comments:

Post a Comment