- 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 hole0
is for it to be in hole1
then move to hole0
, in other wordsP[0] = P[1]*0.5
- If
j = N-1
thenP[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 holej
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