Quantum Grasshopper

Probability ; Difficulty: hard

Can you catch the Grasshopper?

This week’s Riddler Classic features a grasshopper whose position on a log is only probabilistically determined (a quantum grasshopper?). The Classic reproduced:

You are trying to catch a grasshopper on a balance beam that is 1 meter long. Every time you try to catch it, it jumps to a random point along the interval between 20 centimeters left of its current position and 20 centimeters right of its current position.

If the grasshopper is within 20 centimeters of one of the edges, it will not jump off the edge. For example, if it is 10 centimeters from the left edge of the beam, then it will randomly jump to anywhere within 30 centimeters of that edge with equal probability (meaning it will be twice as likely to jump right as it is to jump left).

After many, many failed attempts to catch the grasshopper, where is it most likely to be on the beam? Where is it least likely? And what is the ratio between these respective probabilities?

*spoiler alert* solution and answer follow below

Building intuition

There are more starting positions available to the grasshopper to hit the central region of the log vs. the edge. In fact there are twice as many places it can start from to land on the central region (20 cm away from each edge) vs. an edge. So intuitively, after a long series of jumps we expect the this region to be the most likely location: one where the grasshopper is twice as likely to be found vs the ends, the least likely positions. Let us show that this is indeed the case.

Building the solution:

Our grasshopper’s jumps are at discrete time periods but its jump positions are not discrete. This makes this a Markov Chain process—one that is still discrete time but with a continuous state space.

Consequently, in place of one-step transition probabilities, we now have to deal with analogous one-step transition densities (pdf).

Let \(p(x,y)\) be the transition density function: i.e. given a definite starting position \(x\) in cms measured from one end, the probability density of the grasshopper being at \(y\) cms after its next hop. Given the hop parameters, we can write:

\begin{align*}
p(x,y) &= \frac{1}{x + 20} && 0\leq x \leq 20 &&&\textit{ and } &&&&0\leq y \leq x+20 \\
&= \frac{1}{40} && 20\leq x \leq 80 &&&\textit{ and } &&&& x-20\leq y \leq x+20 \\
&= \frac{1}{100 -x + 20} \quad && 80\leq x \leq 100 &&&\textit{ and } &&&& x- 20\leq y \leq 100 \\
&= 0 &&\textit{everywhere else}
\end{align*}

A regular irreducible and aperiodic Discrete Time Markov Chain is well behaved and has a limiting and stationary distribution: a probability distribution over states that it eventually reaches regardless of the initial/starting conditions.

Our analogous continuous space model also has such a limiting/stationary density distribution \(\pi\) of the grasshopper on the log such that any next hop results in an identical distribution. i.e.

\[
\pi(y) = \int_{S}^{}\pi(x) \ p(x,y) \ dx
\]

where S is the state space which in our case is the length of the log! This is the distribution we seek as an answer to the question.

In the discrete space case, linear algebra offers a simple way to calculate the stationary or limiting distribution. How do we go about finding this stationary distribution function \(\pi(x)\) in our continuous state space?

Let’s see what we can learn by converting the problem to the more familiar discrete space equivalent.

The discrete case:

Let’s divide the log into ten 10cm intervals and restrict the grasshopper’s jumps to one of these 11 points only, with a jump range of at most 2 points or 20 cms as in the continuous case.

This now becomes a discrete time (and space) Markov chain with 11 states and a transition probability matrix:

\[
P = \begin{bmatrix}
0 & 0.5 & 0.5 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0.33 & 0 & 0.33 & 0.33 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0.25 & 0.25 & 0 & 0.25 & 0.25 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0.25 & 0.25 & 0 & 0.25 & 0.25 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0.25 & 0.25 & 0 & 0.25 & 0.25 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0.25 & 0.25 & 0 & 0.25 & 0.25 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0.25 & 0.25 & 0 & 0.25 & 0.25 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0.25 & 0.25 & 0 & 0.25 & 0.25 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0.25 & 0.25 & 0 & 0.25 & 0.25 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0.33 & 0.33 & 0 & 0.33 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0.5 & 0.5 & 0 \\
\end{bmatrix}
\]

Since this is finite, irreducible and aperiodic, the steady state stationary distribution vector \(\pi\) is readily obtained as the solution to:

\[
\pi\cdot P = \pi \\
\sum_{i}^{}\pi_i = 1
\]

Some simple code (shared at the end) gives us the solution shown in the chart.

This discrete case example gives us the needed intuition of what the continuous case stationary density function could look like! Uniform between 20cm and 80cm and declining linearly to half the value towards the ends.

This intuition is further reinforced by trying this out with 10001 points on the log and jumps between 1 and 200 points. This is shown in the chart below!

Final solution and proof!

It’s clear then that the long run distribution to attempt for the continuous case should follow the same trapezoidal shape—uniform between 20cm and 80 cm and decreasing linearly to half this uniform value at the ends.

The constraints above plus the condition that the grasshopper has to be somewhere on the log gives us after some simple math, the following steady state continuous probability density function (p.d.f) for the grasshopper:

\begin{align*}
\pi(x) &= \frac{1}{180} \ + \ \frac{x}{20\cdot 180 } \qquad && 0\leq x \leq 20 \\
&= \frac{1}{90} && 20\leq x \leq 80 \\
&= \frac{1}{180} \ + \ \frac{100-x}{20\cdot 180 } \quad && 80\leq x \leq 100 \\
&= 0 &&\textit{everywhere else – grasshopper has to be on the log!}
\end{align*}

It’s a straightforward exercise to verify that the area under the pdf is 1 as it should for a valid pdf. The probability density (pdf) is shown below.

Clearly the max to the min is 1/90 divided by 1/180 = 2.

But wait, you might say: Where’s the proof that this is indeed a stationary pdf for the continuous state space case? Fair point. This is proved in the next section below the pdf figure.

To prove that the pdf \(\pi (x)\) is a stationary distribution, we have to show:

\[
\pi(y) = \int_{S}^{}\pi(x) \ p(x,y) \ dx
\]

This just says that for a distribution to be steady state, a transition from this state has to result in an identical distribution for the resulting state. That’s what makes it a steady state! Let’s show this for \( 0\leq y\leq 20 \).

For a hop to land between 0 and 20 cms, it has to start between 0 and \(y + 20\). i.e. \(x\) has to be between 0 and 20 OR between 20 and \(y + 20\). Using \(\pi(x)\) from above and \(p(x,y)\) from the section at the beginning:

\begin{align*}
\pi(y) &= \int_{S \equiv 0\leq y \leq 20 }^{}\pi(x) \ p(x,y) \ dx \\
\pi(y) &= \int_{0}^{20} \left [ \frac{1}{180} \ + \ \frac{x}{20\cdot 80}\right ] \frac{1}{x+20} \ dx \ \ + \ \ \int_{20}^{y+20} \frac{1}{90}\cdot \frac{1}{40}\ dx \\
&= \frac{1}{180} + \frac{y}{20\cdot 180} \qquad (0\leq y\leq 20) \\ \\
&= \textit{same as } \pi(x) \textit{ for the range [0,20]}
\end{align*}

The distributions for landing on the other sections of the log can be proved similarly showing that the pdf \(\pi (x)\) is indeed the stationary distribution we sought i.e. this distribution will eventually be reached and a transition from this distribution will result in a state with identical pdf making it a steady state distribution

Final Observations:

As the grasshopper’s jump varies between 1 and 50 cms, the general shape of the trapezoid stays the same; the point at which it switches from being uniform shifts to the new jump points. But the overall max to min probability (center to edge stays the same at 2)
For the grasshopper location to have a uniform probability distribution, the jump range has to be the length of the log = 100 cm as you would expect.
One wonders if this grasshopper can pass through two slits at once 🤔.

Code to calculate steady state probabilities in the discrete case:

				
					import numpy as np 

num_states = 101          #number of discrete state points
jump_range = 20           #jump range 

tp = []                   #tansition probability matrix

for i in range (0,num_states,1):
	if i + jump_range > num_states - 1:
		upper = num_states - 1
	else:
		upper = i + jump_range

	if i - jump_range < 0:
		lower = 0
	else:
		lower = i - jump_range

	p = 1/(upper - lower)	
	tp.append([])

	for j in range (0, num_states,1):
		if abs(j-i) > jump_range or j == i:
			tp[i].append(0)
			
		else:
			tp[i].append(p)


# solving the steady state equation

p_hat = tp.copy()

for i in range(0,num_states,1):
	p_hat[i].pop()

b = p_hat.pop()
 
identity = np.identity(num_states-1)

i_p = identity - p_hat
i_p_inv = np.linalg.inv(i_p)
pi_hat = np.dot(b,i_p_inv)

c = 1 

for i in pi_hat:
	c += i 

pi = []       # this is the vector of final solution probabilities

for i in pi_hat:
	pi.append(i/c)

pi.append(1/c)

# print out steady state state probabilities
for i in pi:
	print(i)