Objectives

  • Expand your R programming skills by creating some increasingly sophisticated random walk procedures!

Walk 1 - Grid-Based Walk

The easiest kind of random walk to program occurs along a rectangular grid. It has the following properties and constraints:

  • Ther walker has no memory of its previous state.
  • Step length is limited to 1 unit.
  • Direction is limited to vertical or horizontal.
  • Each time step, the walker may move a single step in either the horizontal or vertical direction.

Walk Code 1

Here’s some code to perform a simple grid-based random walk:

# How many steps to walk?
n_steps = 400

# Empty matrix to keep track of the position
data_steps = matrix(0, ncol = 2, nrow = n_steps + 1)

# Starting coordinates
x0 = 0; y0 = 0

# before step, set x1 y1 to x0, y0
x_i = x0; y_i = y0


# Loop for the walk:
for(i in 1:n_steps)
{
  direction_i = rbinom(1, 1, 0.5)
  # horizontal
  if (direction_i == 0) {
    x_i = x_i + sample(c(-1, 1), 1)
  }
  # vertical
  if (direction_i == 1)
  {
    y_i = y_i + sample(c(-1, 1), 1)
  }
  data_steps[i+1, 1] = x_i
  data_steps[i+1, 2] = y_i
}

Walk 1 Plot

Now let’s plot it!

par(mar = c(4, 4, 1, 1))
plot(
  data_steps[, 1],
  data_steps[, 2],
  type = "l", xlab = "x", ylab = "y",
  asp = 1)

It would be helpful to show the origin and endpoints:

par(mar = c(4, 4, 1, 1))
plot(
  data_steps[, 1],
  data_steps[, 2],
  type = "l", xlab = "x", ylab = "y",
  asp = 1)

points(
  data_steps[1, 1],
  data_steps[1, 2],
  pch = 16, col = 1, cex = 2)
points(
  data_steps[n_steps + 1, 1],
  data_steps[n_steps + 1, 2],
  pch = 16, col = 2, cex = 2)

Walk 1: Function version

A single random walk is cool, but the true power of random walks is unlocked when we can observe the behavior of multiple walks.

Whenever we want to perform a task more than once, it’s a good idea to wrap our code into a function.

I could easily adapt the code above to fit inside a custom function. I’ll provide the function specification, and you can fill in the details!

r_walk_grid = function(
    n_steps,
    step_length = 1,
    x0 = 0, y0 = 0)
{
  # your code goes here 
}

Your function should return a matrix containing the simulation data.

To test out your function, you can run the following code:

walk_1_dat = r_walk_grid(n_steps = 400)

plot(walk_1_dat[, 1], walk_1_dat[, 2], type = "l", xlab = "x", ylab = "y", asp = 1)
points(
  walk_1_dat[1, 1],
  walk_1_dat[1, 2],
  pch = 16, col = 1, cex = 2)
points(
  walk_1_dat[n_steps + 1, 1],
  walk_1_dat[n_steps + 1, 2],
  pch = 16, col = 2, cex = 2)

Walk 1 - Multiple Realizations

Now that you have a function, it’s relatively easy to create multiple simulations and plot them!

Here’s a very simple loop-based method you might use to create and plot several realizations

n_steps = 400
n_sims = 20

walk_1_dat = r_walk_grid(n_steps = n_steps)

plot(
  walk_1_dat[, 1],
  walk_1_dat[, 2],
  type = "l", asp = 1,
  xlab = "x", ylab = "y",
  xlim = c(-50, 50),
  ylim = c(-50, 50))
points(
  walk_1_dat[1, 1],
  walk_1_dat[1, 2],
  pch = 16, col = 1, cex = 2)
points(
  walk_1_dat[n_steps + 1, 1],
  walk_1_dat[n_steps + 1, 2],
  pch = 16, col = 2, cex = 2)

for (i in 1:(n_sims - 1))
{
  walk_1_dat = r_walk_grid(n_steps = n_steps)
  
  points(
    walk_1_dat[, 1],
    walk_1_dat[, 2],
    type = "l")
  points(
    walk_1_dat[1, 1],
    walk_1_dat[1, 2],
    pch = 16, col = 1, cex = 2)
  points(
    walk_1_dat[n_steps + 1, 1],
    walk_1_dat[n_steps + 1, 2],
    pch = 16, col = 2, cex = 2)
  
}

Plot Function


Our code to plot the random walk, as well as the code to run the loop are already quite complicated. You could easily create additional functions to run the simulation loop and do the plotting.

For example, I can write a simple function to plot a simulation run:

plot_sim = function(
    dat_sim,
    add = FALSE,
    ylab = "y", 
    xlab = "x",
    col_start = 1, col_end = 2,
    cex_start = 2, cex_end = 2,
    pch_start = 16, pch_end = 16,
    ...)
{
  
  # If this is the first simulation, create a new plot
  if (!add)
  {
    plot(
      dat_sim[, 1], dat_sim[, 2],
      xlab = xlab, ylab = ylab,
      type = "l", asp = 1,
      ...)
  }
  
  # If adding to an existing plot:
  if (add)
  {
    lines(
      dat_sim[, 1], dat_sim[, 2],
      ...)
  }
  
  points(
    dat_sim[1, 1], dat_sim[1, 2],
    pch = pch_start, col = col_start,
    cex = cex_start)
  points(
    dat_sim[nrow(dat_sim), 1],
    dat_sim[nrow(dat_sim), 2],
    pch = pch_end, col = col_end,
    cex = cex_end)
}

Now to test it out:

n_steps = 500
n_sims = 10

xlm = c(-50, 50)
ylm = c(-50, 50)

walk_1_dat = r_walk_grid(n_steps = n_steps)

plot_sim(walk_1_dat, xlim = xlm, ylim = ylm, lwd = 2)


for (i in 1:n_sims)
{
  dat_i = r_walk_grid(n_steps = n_steps)
  plot_sim(dat_i, xlim = xlm, ylim = ylm, lwd = 2, add = TRUE)
}

It could use some improvements, for example the lines for later simulations cover up the points for early sims, but it’ll do for now. Feel free to add your own improvements!

Walk 2 - Continuous Direction

For the first elaboration, let’s allow our walker to turn in any direction, i.e. not be restricted to walking on a grid.

Trigonometry Review

We’ll need some verfy simple trigonometry for this part (sorry).

The important things to recall are:

  • the sin and cos functions require a single argument: the angle.
  • cos() gives us the x-coordinate.
  • sin() gives us the y-coordinate.
  • The trig functions in R expect angles to be in radians.
  • There are \(2\pi\) radians in a complete circle.
  • This means that the angles \(3\pi\) and \(\pi\) are equivalent. (why?)

Walk Code 2

For this random walk implementation, let’s keep the step length constant, but allow the step direction to vary uniformly over the range \(0\) to \(2\pi\).

Fortunately our code only needs some minor modifications to accommodate the continuous direction enhancement.

I’ll let you create your own implementation.

  • Your implementation should via a function called r_walk_angle().
  • It should take the same arguments as r_walk_grid()

Your

Test your function

You can compare your function to my reference implementation graphically:

plot_sim(r_walk_angle(5))

Things to note:

  • The constant step length (this is harder to see with lots of steps and multiple simulations ont he same plot).
  • The random direction.

Here’s a realization of multiple simulations:

plot_sim(r_walk_angle(5), xlim = c(-5, 5), ylim = c(-5, 5))
for (i in 1:10)
  plot_sim(r_walk_angle(5), xlim = c(-5, 5), ylim = c(-5, 5), add = T)

And another with lots of simulations:

plot_sim(r_walk_angle(5), xlim = c(-5, 5), ylim = c(-5, 5))
for (i in 1:1000)
  plot_sim(
    r_walk_angle(5),
    xlim = c(-5, 5),
    ylim = c(-5, 5),
    add = T, 
    lwd = 0.1, 
    col_end = adjustcolor("red", 0.2))

Walk 3: Uniformly-Distributed Step Lengths

For this iteration, let’s allow the step length to vary according to a uniform distribution.

This is another easy modification that you can implement in a function.

  • You’ll need some new arguments to your function: minimum and maximum step lengths for the random uniform numbers.

The function spec should look like this:

r_walk_angle_unif = function(
    n_steps,
    step_min = 0.1,
    step_max = 1,
    x0 = 0, y0 = 0)
{
  # Your implementation here.
}

The results of my implementation looked like:

plot_sim(r_walk_angle_unif(10))

Walk 4: Correlated Random Walk

For the final walk, you’ll implmement a correlated random walk.

You might ask: what’s a correlated random walk?

A correlated random walk is a walk in which the walker remembers the direction of its most recent step. The direction of the next step is constrained by the previous direction. For example, a walker may be limited to a maximum turn of 90 degrees (\(\pi / 4\) radians).

For a correlated random walk, you have a number of choices to make, including:

  • How wide should the ‘cone’ of allowed turing angles be?
  • Should the turning angle be uniform within the cone?

For my implementation, I chose to allow the turning angle to vary uniformly within the cone.

The spec for my function was:

r_walk_angle_correlated = function(
    n_steps,
    step_min = 0.1,
    step_max = 1,
    cone_width = pi/4,
    initial_direction = 0,
    x0 = 0, y0 = 0)
{
  # Your implementation here 
}

Here’s how my simulation looked with 1000 steps:

plot_sim(r_walk_angle_correlated(1000))

Here’s how my simulation looked with 1000 steps and multiple realizations:

plot_sim(r_walk_angle_correlated(1000), xlim = c(-400, 400), ylim = c(-400, 400))

for (i in 1:15)
  plot_sim(
    r_walk_angle_correlated(1000),
    xlim = c(-5, 5),
    ylim = c(-5, 5),
    add = T, 
    lwd = 0.1, 
    col_end = adjustcolor("red", 0.5))

Note that I could specify a cone width of 360 degrees to make my simulation an uncorrelated random walk:

plot_sim(r_walk_angle_correlated(1000, cone_width = 2 * pi))