Modeling Levy Walk in R

Question

I am trying to create a series of numbers to simulate a Levy Walk in R. I am currently using the following code:

alpha=2 n=1000 x=rep(0,n) y=rep(0,n) for (i in 2:n){ theta=runif(1)*2*pi f=runif(1)^(-1/alpha) x[i]=x[i-1]+f*cos(theta) y[i]=y[i-1]+f*sin(theta) } 

The code works as expected, and I can generate numbers according to my requirements. The figure below shows such a “Levy Hall”:

The following histogram confirms that the numbers generated (i.e. f) actually belong to the power law:

My question is this: The resulting step lengths (i.e. F) are quite large. Can I change the code so that the step lengths fall only at some boundary [fmin, fmax]?

PS I intentionally did not vectorize the code.

Answer 1

Try using this:

 f=runif(1, fmax^(-alpha), fmin^(-alpha))^(-1/alpha) 

Note that you need 0 < fmin < fmax .

By the way, you can vectorize your code as follows:

 theta <- runif(n-1)*2*pi f <- runif(n-1, fmax^(-alpha), fmin^(-alpha))^(-1/alpha) x <- c(0, cumsum(f*cos(theta))) y <- c(0, cumsum(f*sin(theta))) 

Answer 2

Just for accuracy, what you pretend to be here is a Levy flight. For this to be a walk in Levy, you must allow the particle to “walk” from the beginning to the end of each flight (for example, using for ). If you build the final simulation using plot(x, y, type = "o") , you will see that with the code there are no positions in flight (without walking).