numpy
has useful functions for generating random numbers and maninulating arrays.
matplotlib
has useful functions for plotting things.
import numpy as np
import matplotlib.pyplot as plt
We will explicitly model the process of picking a bridge and then picking a creature, rejecting any trial in which a gnome is picked.
numpy.random.choice
is a convenient way to pick from a set of objects with specified weights (i.e draw from a categorical distribution).
bridges = ['TGGGG', 'TTGGG', 'GGGGG']
ps = [0.2, 0.2, 0.6]
trials = []
while len(trials) < 100000:
bridge = np.random.choice(bridges, p=ps)
creature = np.random.choice(list(bridge))
if creature == 'T':
if bridge.count('T') == 1:
trials.append(1)
else:
trials.append(0)
trials
is now a list of 0's and 1's that we need to convert into a list of "fraction of safe passages so far".
cumulative_fractions = np.cumsum(trials) / (np.arange(len(trials)) + 1.)
plt.semilogx(cumulative_fractions)
plt.ylim(0, 1)
plt.axhline(1. / 3, color='red', linestyle='--')
plt.xlabel('Number of trials')
plt.ylabel('Cumulative fraction of safe passages');
The simulation takes a while, in part because we need to generate about $1 / (0.2 * 0.2 + 0.4 * 0.2) \approx 8.3$ trials per troll captured (can you see why?) and in part because np.random.choice
is just slow.
%timeit np.random.choice(bridges, p=ps)
10000 loops, best of 3: 143 µs per loop