Roberto De Almeida

This is a solution to this problem:

A drunken sailor returns to his ship via a plank 15 paces long and 7 paces wide. With each step he has an equal chance of stepping forward, left, right, or standing still. What is the probability that he returns safely to his ship?

The article shows two different solutions, one using a Monte Carlo approach and the other, Markov chains. The Monte Carlo approach is pretty trivial; here’s how I did it using Python:

def walk():
    pos = array([4, 0])
    while 1:
        if pos[1] == 16:
            return True
        elif not 0 < pos[0] < 8:
            return False

        pos += choice([
                array([ 0, 1]),  # forward
                array([ 1, 0]),  # right
                array([-1, 0]),  # left
                array([ 0, 0]),  # still
                ])

This is our random walk. We keep updating the position of the sailor until he either reaches the ship of falls into water. We can estimate the chance that the sailor will reach the boat by running this function several times and counting the percentage that returns True:

REALIZATIONS = 10000
successes = 0
for i in range(REALIZATIONS):
    if walk(): successes += 1
print float(successes)/REALIZATIONS

For 10 thousand realizations I get a probability of 12% that the sailor can cross the plank. But I’m more interested in the Markov chain solution. The idea behind this solution is to create a matrix of the probabilities of going from a place A to a point B. In order to do this I create an array of shape (7+2, 15+1, 7+2, 15+1). This array includes two additional side points, for water, and one additional point along the plank, representing the boat. In this array, the value at the position (i0, j0, i1, j1) represents the probability of a person at point (i0, j0) going to the position (i1, j1). We start with an empty arrays of zeros:

WIDTH = 7
LENGTH = 15
markov = zeros((WIDTH+2, LENGTH+1, WIDTH+2, LENGTH+1), 'f')

Now let’s set the probabilities for the plank:

for i in range(1, WIDTH+1):
    for j in range(LENGTH):
        markov[i, j, i, j+1] = 0.25  # forward
        markov[i, j, i+1, j] = 0.25  # right
        markov[i, j, i-1, j] = 0.25  # left
        markov[i, j, i, j]   = 0.25  # still

So, if the sailor is in the point (i, j) inside the plank it has a 25% chance of moving to point (i+1, j), ie, to the right. We also set the probabilities for when the sailor is on the water; in this case, he will stay there:

# when in water, stay
for j in range(LENGTH+1):
    markov[0, j, 0, j] = 1
    markov[WIDTH+1, j, WIDTH+1, j] = 1

And we do the same when he reaches the boat:

# when in the ship, stay
for i in range(WIDTH+2):
    markov[i, LENGTH, i, LENGTH] = 1

We also reshape our array, making it 2D. It will still keep the probabilistic correspondence between point A and B, but now our space dimension has been turned into 1D:

markov.shape = ((WIDTH+2)*(LENGTH+1), (WIDTH+2)*(LENGTH+1))

The rest now is easy. The transition matrix we just built has the nice property that it represents the chance of going from one point to the other, and if we keep multiplying it by itself n times we have the chance of reaching point B after n steps from point A. So we just need to convert our array into a matrix (there’s a difference in Numpy) and take the power to a resonable number of steps:

MAXSTEPS = 200
result = matrix(markov)**MAXSTEPS
result = array(result).reshape(WIDTH+2, LENGTH+1, WIDTH+2, LENGTH+1)

What we want is the chance of a sailor leaving from point (4, 0) to reach the boat at any point in (i, 15). We can get this easily:

print sum(result[4, 0, :, 15])

And the result is 15%. Pretty cool.