- There are 100 holes in a line, and there is a bunny hiding in one of the holes
- You can only look in one hole at a time, and every time you look the bunny randomly jumps to an adjacent hole

The *intended* solution is some boring algorithm I'm not going to mention here, we're going to explore how to approach this from a probability theory standpoint (similar to Solving Wordle using information theory except my approach doesn't work)

We're going to keep a probability distribution over where the bunny could be, we start off with equal likelyhood for each hole

```
N = 100
P = [1/N] * N
```

Now we write the function to update our belief when we learn the bunny isn't in hole `j`

```
def guess(P, j):
# after setting P[j] = 0 we have removed P[j] probability mass from
# the distribution, so the new sum is 1-P[j]. Dividing everything by
# 1-P[j] fixes this. (a minor detail, instead of 1 we use sum(P) for
# numerical stability since sum(P) isn't exactly 1)
s = sum(P) - P[j]
return [
0 if i == j else P[i] / s
for i in range(N)
]
```

Now for the tricky bit, what happens to our beliefs when the bunny jumps to an adjacent hole?

It might seem hopelessly complicated, but you can figure it out by conditioning on three different cases

- If
`j = 0`

then the only way for the bunny to be in hole`0`

is for it to be in hole`1`

then move to hole`0`

, in other words`P[0] = P[1]*0.5`

- If
`j = N-1`

then`P[N-1] = P[N-2]*0.5`

by the same logic - Otherwise
`P[j] = P[j-1]*0.5 + P[j+1]*0.5`

since we could get to hole`j`

from the left or the right.

This is basically code, the one thing to be careful of is mutation- If we modified `P[j]`

implace it would screw up the computation of `P[j+1]`

. Hence we make a copy

```
def jumps(P):
Pn = P.copy()
Pn[0] = P[1]*0.5
Pn[N-1] = P[N-2]*0.5
for j in range(1, N-1):
Pn[j] = 0.5*P[j-1] + 0.5*P[j+1]
return Pn
```

We're pretty much done, we can wrap it up to test it

```
import random
def argmax(l):
return max(range(len(l)), key=lambda i: l[i])
hole = random.randint(0, N-1)
t = 0
while True:
t += 1
j = argmax(P)
if j == hole:
print(f'Found the bunny in hole {j} after {t} tries!')
break
else:
print(f'[{t}] Guessed {j} with P[{j}] = {P[j]:.4f} (true is {hole})')
# wrong hole, update our beliefs
P = guess(P, j)
P = jumps(P)
# the bunny moves!
if hole == 0: hole += 1
elif hole == N-1: hole -= 1
else: hole += random.choice([-1, 1])
```

Now we can test it

```
$ python3 bunny.py
[1] Guessed 0 with P[0] = 0.0100 (true is 93)
[2] Guessed 2 with P[2] = 0.0101 (true is 92)
[3] Guessed 4 with P[4] = 0.0102 (true is 93)
...
[42] Guessed 82 with P[82] = 0.0176 (true is 86)
[43] Guessed 84 with P[84] = 0.0179 (true is 87)
Found the bunny in hole 86 after 44 tries!
```

Cool, I wonder how well it does on average? If we factor the above loop into a function we can do something like

```
r = [play() for _ in range(10000)]
import statistics as stats
mean, median, stdev = stats.mean(r), stats.median(r), stats.stdev(r)
print(f'mean: {mean} median: {median} stddev {stdev}')
```

Running this I get

`mean: 94.7903 median: 53.0 stddev 118.2896460809102`

How does this compare to random guessing? replacing the `j = argmax(P)`

line with `j = random.randint(0, N-1)`

gives

`mean: 99.5563 median: 70.0 stddev 97.57317399130908`

Our algorithm barely beats random guessing :( something must have gone wrong, I'll figure it out later

The mean for random guessing being close to `100`

makes sense since every time we guess randomly we have a `p = 1/100`

(independent) chance of success, which is exactly the Geometric Distribution meaning the expected value is `1/p = 100`

lining up with our numerical result of around `99.6`

.

PS: Full code is downloadable here