Skip to content

Instantly share code, notes, and snippets.

@pixyj
Created May 7, 2015 14:08
Show Gist options
  • Save pixyj/e1da1568abd2522fae82 to your computer and use it in GitHub Desktop.
Save pixyj/e1da1568abd2522fae82 to your computer and use it in GitHub Desktop.
A solution to a probability problem using a Monte Carlo simulation
from __future__ import division
import random
import itertools
"""
People are waiting in line to board a 100-seat airplane. Steve is the first person in the line. He gets on the plane but suddenly can't remember what his seat number is, so he picks a seat at random. After that, each person who gets on the plane sits in their assigned seat if it's available, otherwise they will choose an open seat at random to sit in.
The flight is full and you are last in line. What is the probability that you get to sit in your assigned seat?
Usage:
python monte.py
"""
def run(people, p, trials):
trial_prob = []
for i in range(trials):
tp = trial(people, p)
trial_prob.append(tp)
#return trial_prob
return sum(trial_prob) / len(trial_prob)
def trial(people, p):
status = {}
for i in range(people):
status[i] = False
choices = [(True, p), (False, 1 - p)]
result = weighted_choice(choices)
if result:
#print "Steve remembered. So he chooses 0"
status[0] = True
else:
position = choose_random_position_among_available(status)
status[position] = True
#print "Steve forgot and he chooses {}".format(position)
for i in range(1, people-1):
if not status[i]:
status[i] = True
#print "Person {} enters the room. Her seat is NOT occupied. She chooses {}".format(i, i)
else:
position = choose_random_position_among_available(status)
status[position] = True
#print "Person {} enters the room. Her set is OCCUPIED. So she chooses {}".format(i, position)
return 1 if status[people-1] == False else 0
def get_not_occupied_positions(status):
positions = []
for position, value in status.iteritems():
if value == False:
positions.append(position)
return positions
def choose_random_position_among_available(status):
not_occupied_positions = get_not_occupied_positions(status)
position = random.choice(not_occupied_positions)
return position
def weighted_choice(choices):
total = sum(w for c, w in choices)
r = random.uniform(0, total)
upto = 0
for c, w in choices:
if upto + w > r:
return c
upto += w
assert False, "Shouldn't get here"
def run_monte_and_store_results():
results = []
for i in range(0,101):
p = i / 100
result = run(100, p, 10000)
results.append("{},{}".format(p, result))
print "Progress: {}%".format(i)
s = '\n'.join(results)
with open("monte.csv", "w") as f:
f.write(s)
if __name__ == '__main__':
run_monte_and_store_results()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment