%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import numpy as np
import pandas as pd
from scipy.stats import beta, gamma, binom
import matplotlib.pyplot as plt
plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = 9, 6
A clinic has three doctors. Patients come into the clinic at random, starting at 9 a.m., according to a Poisson process with time parameter 10 minutes: that is, the time after opening at which the first patient appears follows an exponential distribution with expectation 10 minutes and then, after each patient arrives, the waiting time until the next patient is independently exponentially distributed, also with expectation 10 minutes. When a patient arrives, he or she waits until a doctor is available. The amount of time spent by each doctor with each patient is a random variable, uniformly distributed between 5 and 25 minutes. The office stops admitting new patients at 4 p.m. and closes when the last patient is through with the doctor.
Simulate this process once. How many patients came to the office? How many had to wait for a doctor? What was the average wait time of those who did have to wait? What was the longest wait? How many minutes past 4 p.m. did the office close?
Simulate the process 2000 times and estimate the quartiles for each of the sum- maries in part (a). Summarize your results in a 5 × 3 table with the median in the middle column.
Use your 2000 simulations to estimate the probability, on a particular day, that
n_minutes = 420
class Patient:
def __init__(self):
self.wait_time = 0
self.waited = False
def __str__(self):
return '<patient: wait={}>'.format(self.wait_time)
def __repr__(self):
return self.__str__()
class Doctor:
def __init__(self):
self.time_left = None
self.is_free = True
def take_appt(self):
self.time_left = np.ceil(np.random.uniform(5,25)).astype(int)
self.is_free = False
def __str__(self):
if self.is_free:
return '<Doc: free>'
return '<Doc: time left={}>'.format(self.time_left)
def __repr__(self):
return self.__str__()
def __nonzero__(self):
return self.is_free
def __bool__(self):
return self.is_free
def get_arrival_times(rate=10):
times = []
t = int(np.random.exponential(rate))
while t < n_minutes:
times.append(t)
t = t + int(np.random.exponential(rate))
return frozenset(times)
def run_simulation(time_param=10, verbose=False):
def get_open_doc():
for doc in doctors:
if doc.is_free:
return doc
arrival_times = get_arrival_times(time_param)
doctors = [Doctor() for _ in range(3)]
patients = []
wait_times = []
stats = {
'n_patients': 0,
'n_waiters': 0
}
# run the simulation
t = 0
while t < n_minutes or not all(doctors): # in case any doctors are working
if verbose:
print(t)
# add a patient upon arrival
if t in arrival_times and t < n_minutes:
if verbose:
print('patient arrived!')
stats['n_patients'] += 1
patients.append(Patient())
# decrement counters, free up doctors whom are done
for doc in doctors:
if not doc.is_free:
doc.time_left -= 1
if not doc.is_free and doc.time_left <= 0:
doc.is_free = True
doc.time_left = None
# send patient to doctor if one is open
while any(doctors) and patients:
p = patients.pop(0)
wait_times.append(p.wait_time)
if p.waited:
stats['n_waiters'] += 1
d = get_open_doc()
d.take_appt()
# increment wait times
for patient in patients:
patient.wait_time += 1
if not patient.waited:
patient.waited = True
if verbose:
print('patients:', patients)
print('doctors: ', doctors)
print()
t += 1
stats['avg_wait'] = np.mean([x for x in wait_times if x != 0]) if any(wait_times) else 0
stats['longest_wait'] = max(wait_times)
stats['mins_past_close'] = t - 420
stats['wait_times'] = wait_times
return stats
run_simulation(verbose=True)
0 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 1 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 2 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 3 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 4 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 5 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 6 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 7 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 8 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 9 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 10 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 11 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 12 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 13 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 14 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 15 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 16 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 17 patient arrived! patients: [] doctors: [<Doc: time left=22>, <Doc: free>, <Doc: free>] 18 patients: [] doctors: [<Doc: time left=21>, <Doc: free>, <Doc: free>] 19 patients: [] doctors: [<Doc: time left=20>, <Doc: free>, <Doc: free>] 20 patients: [] doctors: [<Doc: time left=19>, <Doc: free>, <Doc: free>] 21 patients: [] doctors: [<Doc: time left=18>, <Doc: free>, <Doc: free>] 22 patients: [] doctors: [<Doc: time left=17>, <Doc: free>, <Doc: free>] 23 patients: [] doctors: [<Doc: time left=16>, <Doc: free>, <Doc: free>] 24 patients: [] doctors: [<Doc: time left=15>, <Doc: free>, <Doc: free>] 25 patients: [] doctors: [<Doc: time left=14>, <Doc: free>, <Doc: free>] 26 patients: [] doctors: [<Doc: time left=13>, <Doc: free>, <Doc: free>] 27 patients: [] doctors: [<Doc: time left=12>, <Doc: free>, <Doc: free>] 28 patients: [] doctors: [<Doc: time left=11>, <Doc: free>, <Doc: free>] 29 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 30 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 31 patients: [] doctors: [<Doc: time left=8>, <Doc: free>, <Doc: free>] 32 patients: [] doctors: [<Doc: time left=7>, <Doc: free>, <Doc: free>] 33 patients: [] doctors: [<Doc: time left=6>, <Doc: free>, <Doc: free>] 34 patients: [] doctors: [<Doc: time left=5>, <Doc: free>, <Doc: free>] 35 patient arrived! patients: [] doctors: [<Doc: time left=4>, <Doc: time left=7>, <Doc: free>] 36 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=6>, <Doc: free>] 37 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=5>, <Doc: free>] 38 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=4>, <Doc: free>] 39 patient arrived! patients: [] doctors: [<Doc: time left=22>, <Doc: time left=3>, <Doc: free>] 40 patient arrived! patients: [] doctors: [<Doc: time left=21>, <Doc: time left=2>, <Doc: time left=19>] 41 patient arrived! patients: [<patient: wait=1>] doctors: [<Doc: time left=20>, <Doc: time left=1>, <Doc: time left=18>] 42 patients: [] doctors: [<Doc: time left=19>, <Doc: time left=24>, <Doc: time left=17>] 43 patients: [] doctors: [<Doc: time left=18>, <Doc: time left=23>, <Doc: time left=16>] 44 patients: [] doctors: [<Doc: time left=17>, <Doc: time left=22>, <Doc: time left=15>] 45 patients: [] doctors: [<Doc: time left=16>, <Doc: time left=21>, <Doc: time left=14>] 46 patients: [] doctors: [<Doc: time left=15>, <Doc: time left=20>, <Doc: time left=13>] 47 patients: [] doctors: [<Doc: time left=14>, <Doc: time left=19>, <Doc: time left=12>] 48 patients: [] doctors: [<Doc: time left=13>, <Doc: time left=18>, <Doc: time left=11>] 49 patients: [] doctors: [<Doc: time left=12>, <Doc: time left=17>, <Doc: time left=10>] 50 patients: [] doctors: [<Doc: time left=11>, <Doc: time left=16>, <Doc: time left=9>] 51 patients: [] doctors: [<Doc: time left=10>, <Doc: time left=15>, <Doc: time left=8>] 52 patients: [] doctors: [<Doc: time left=9>, <Doc: time left=14>, <Doc: time left=7>] 53 patients: [] doctors: [<Doc: time left=8>, <Doc: time left=13>, <Doc: time left=6>] 54 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=12>, <Doc: time left=5>] 55 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=11>, <Doc: time left=4>] 56 patients: [] doctors: [<Doc: time left=5>, <Doc: time left=10>, <Doc: time left=3>] 57 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=9>, <Doc: time left=2>] 58 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=8>, <Doc: time left=1>] 59 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=7>, <Doc: free>] 60 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=6>, <Doc: free>] 61 patients: [] doctors: [<Doc: free>, <Doc: time left=5>, <Doc: free>] 62 patient arrived! patients: [] doctors: [<Doc: time left=21>, <Doc: time left=4>, <Doc: free>] 63 patients: [] doctors: [<Doc: time left=20>, <Doc: time left=3>, <Doc: free>] 64 patients: [] doctors: [<Doc: time left=19>, <Doc: time left=2>, <Doc: free>] 65 patients: [] doctors: [<Doc: time left=18>, <Doc: time left=1>, <Doc: free>] 66 patients: [] doctors: [<Doc: time left=17>, <Doc: free>, <Doc: free>] 67 patients: [] doctors: [<Doc: time left=16>, <Doc: free>, <Doc: free>] 68 patients: [] doctors: [<Doc: time left=15>, <Doc: free>, <Doc: free>] 69 patients: [] doctors: [<Doc: time left=14>, <Doc: free>, <Doc: free>] 70 patients: [] doctors: [<Doc: time left=13>, <Doc: free>, <Doc: free>] 71 patients: [] doctors: [<Doc: time left=12>, <Doc: free>, <Doc: free>] 72 patients: [] doctors: [<Doc: time left=11>, <Doc: free>, <Doc: free>] 73 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 74 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 75 patient arrived! patients: [] doctors: [<Doc: time left=8>, <Doc: time left=21>, <Doc: free>] 76 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=20>, <Doc: free>] 77 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=19>, <Doc: free>] 78 patients: [] doctors: [<Doc: time left=5>, <Doc: time left=18>, <Doc: free>] 79 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=17>, <Doc: free>] 80 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=16>, <Doc: free>] 81 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=15>, <Doc: free>] 82 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=14>, <Doc: free>] 83 patient arrived! patients: [] doctors: [<Doc: time left=11>, <Doc: time left=13>, <Doc: free>] 84 patients: [] doctors: [<Doc: time left=10>, <Doc: time left=12>, <Doc: free>] 85 patients: [] doctors: [<Doc: time left=9>, <Doc: time left=11>, <Doc: free>] 86 patients: [] doctors: [<Doc: time left=8>, <Doc: time left=10>, <Doc: free>] 87 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=9>, <Doc: free>] 88 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=8>, <Doc: free>] 89 patient arrived! patients: [] doctors: [<Doc: time left=5>, <Doc: time left=7>, <Doc: time left=15>] 90 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=6>, <Doc: time left=14>] 91 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=5>, <Doc: time left=13>] 92 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=4>, <Doc: time left=12>] 93 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=3>, <Doc: time left=11>] 94 patients: [] doctors: [<Doc: free>, <Doc: time left=2>, <Doc: time left=10>] 95 patients: [] doctors: [<Doc: free>, <Doc: time left=1>, <Doc: time left=9>] 96 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=8>] 97 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=7>] 98 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=6>] 99 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=5>] 100 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=4>] 101 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=3>] 102 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=2>] 103 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=1>] 104 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 105 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 106 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 107 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 108 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 109 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 110 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 111 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 112 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 113 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 114 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 115 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 116 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 117 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 118 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 119 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 120 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 121 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 122 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 123 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 124 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 125 patient arrived! patients: [] doctors: [<Doc: time left=8>, <Doc: free>, <Doc: free>] 126 patients: [] doctors: [<Doc: time left=7>, <Doc: free>, <Doc: free>] 127 patients: [] doctors: [<Doc: time left=6>, <Doc: free>, <Doc: free>] 128 patients: [] doctors: [<Doc: time left=5>, <Doc: free>, <Doc: free>] 129 patients: [] doctors: [<Doc: time left=4>, <Doc: free>, <Doc: free>] 130 patients: [] doctors: [<Doc: time left=3>, <Doc: free>, <Doc: free>] 131 patients: [] doctors: [<Doc: time left=2>, <Doc: free>, <Doc: free>] 132 patients: [] doctors: [<Doc: time left=1>, <Doc: free>, <Doc: free>] 133 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 134 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 135 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 136 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 137 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 138 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 139 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 140 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 141 patient arrived! patients: [] doctors: [<Doc: time left=13>, <Doc: free>, <Doc: free>] 142 patients: [] doctors: [<Doc: time left=12>, <Doc: free>, <Doc: free>] 143 patients: [] doctors: [<Doc: time left=11>, <Doc: free>, <Doc: free>] 144 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 145 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 146 patients: [] doctors: [<Doc: time left=8>, <Doc: free>, <Doc: free>] 147 patients: [] doctors: [<Doc: time left=7>, <Doc: free>, <Doc: free>] 148 patients: [] doctors: [<Doc: time left=6>, <Doc: free>, <Doc: free>] 149 patient arrived! patients: [] doctors: [<Doc: time left=5>, <Doc: time left=13>, <Doc: free>] 150 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=12>, <Doc: free>] 151 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=11>, <Doc: free>] 152 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=10>, <Doc: free>] 153 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=9>, <Doc: free>] 154 patients: [] doctors: [<Doc: free>, <Doc: time left=8>, <Doc: free>] 155 patients: [] doctors: [<Doc: free>, <Doc: time left=7>, <Doc: free>] 156 patients: [] doctors: [<Doc: free>, <Doc: time left=6>, <Doc: free>] 157 patients: [] doctors: [<Doc: free>, <Doc: time left=5>, <Doc: free>] 158 patients: [] doctors: [<Doc: free>, <Doc: time left=4>, <Doc: free>] 159 patients: [] doctors: [<Doc: free>, <Doc: time left=3>, <Doc: free>] 160 patient arrived! patients: [] doctors: [<Doc: time left=19>, <Doc: time left=2>, <Doc: free>] 161 patients: [] doctors: [<Doc: time left=18>, <Doc: time left=1>, <Doc: free>] 162 patients: [] doctors: [<Doc: time left=17>, <Doc: free>, <Doc: free>] 163 patient arrived! patients: [] doctors: [<Doc: time left=16>, <Doc: time left=6>, <Doc: free>] 164 patients: [] doctors: [<Doc: time left=15>, <Doc: time left=5>, <Doc: free>] 165 patients: [] doctors: [<Doc: time left=14>, <Doc: time left=4>, <Doc: free>] 166 patients: [] doctors: [<Doc: time left=13>, <Doc: time left=3>, <Doc: free>] 167 patients: [] doctors: [<Doc: time left=12>, <Doc: time left=2>, <Doc: free>] 168 patients: [] doctors: [<Doc: time left=11>, <Doc: time left=1>, <Doc: free>] 169 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 170 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 171 patient arrived! patients: [] doctors: [<Doc: time left=8>, <Doc: time left=7>, <Doc: free>] 172 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=6>, <Doc: free>] 173 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=5>, <Doc: free>] 174 patients: [] doctors: [<Doc: time left=5>, <Doc: time left=4>, <Doc: free>] 175 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=3>, <Doc: free>] 176 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=2>, <Doc: free>] 177 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=1>, <Doc: free>] 178 patients: [] doctors: [<Doc: time left=1>, <Doc: free>, <Doc: free>] 179 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 180 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 181 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 182 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 183 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 184 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 185 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 186 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 187 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 188 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 189 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 190 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 191 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 192 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 193 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 194 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 195 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 196 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 197 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 198 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 199 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 200 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 201 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 202 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 203 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 204 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 205 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 206 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 207 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 208 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 209 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 210 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 211 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 212 patient arrived! patients: [] doctors: [<Doc: time left=22>, <Doc: free>, <Doc: free>] 213 patients: [] doctors: [<Doc: time left=21>, <Doc: free>, <Doc: free>] 214 patients: [] doctors: [<Doc: time left=20>, <Doc: free>, <Doc: free>] 215 patients: [] doctors: [<Doc: time left=19>, <Doc: free>, <Doc: free>] 216 patients: [] doctors: [<Doc: time left=18>, <Doc: free>, <Doc: free>] 217 patients: [] doctors: [<Doc: time left=17>, <Doc: free>, <Doc: free>] 218 patients: [] doctors: [<Doc: time left=16>, <Doc: free>, <Doc: free>] 219 patients: [] doctors: [<Doc: time left=15>, <Doc: free>, <Doc: free>] 220 patients: [] doctors: [<Doc: time left=14>, <Doc: free>, <Doc: free>] 221 patients: [] doctors: [<Doc: time left=13>, <Doc: free>, <Doc: free>] 222 patients: [] doctors: [<Doc: time left=12>, <Doc: free>, <Doc: free>] 223 patients: [] doctors: [<Doc: time left=11>, <Doc: free>, <Doc: free>] 224 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 225 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 226 patients: [] doctors: [<Doc: time left=8>, <Doc: free>, <Doc: free>] 227 patients: [] doctors: [<Doc: time left=7>, <Doc: free>, <Doc: free>] 228 patients: [] doctors: [<Doc: time left=6>, <Doc: free>, <Doc: free>] 229 patients: [] doctors: [<Doc: time left=5>, <Doc: free>, <Doc: free>] 230 patients: [] doctors: [<Doc: time left=4>, <Doc: free>, <Doc: free>] 231 patients: [] doctors: [<Doc: time left=3>, <Doc: free>, <Doc: free>] 232 patients: [] doctors: [<Doc: time left=2>, <Doc: free>, <Doc: free>] 233 patients: [] doctors: [<Doc: time left=1>, <Doc: free>, <Doc: free>] 234 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 235 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 236 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 237 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 238 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 239 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 240 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 241 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 242 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 243 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 244 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 245 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 246 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 247 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 248 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 249 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 250 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 251 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 252 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 253 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 254 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 255 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 256 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 257 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 258 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 259 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 260 patient arrived! patients: [] doctors: [<Doc: time left=25>, <Doc: free>, <Doc: free>] 261 patients: [] doctors: [<Doc: time left=24>, <Doc: free>, <Doc: free>] 262 patients: [] doctors: [<Doc: time left=23>, <Doc: free>, <Doc: free>] 263 patient arrived! patients: [] doctors: [<Doc: time left=22>, <Doc: time left=8>, <Doc: free>] 264 patients: [] doctors: [<Doc: time left=21>, <Doc: time left=7>, <Doc: free>] 265 patients: [] doctors: [<Doc: time left=20>, <Doc: time left=6>, <Doc: free>] 266 patients: [] doctors: [<Doc: time left=19>, <Doc: time left=5>, <Doc: free>] 267 patients: [] doctors: [<Doc: time left=18>, <Doc: time left=4>, <Doc: free>] 268 patients: [] doctors: [<Doc: time left=17>, <Doc: time left=3>, <Doc: free>] 269 patients: [] doctors: [<Doc: time left=16>, <Doc: time left=2>, <Doc: free>] 270 patients: [] doctors: [<Doc: time left=15>, <Doc: time left=1>, <Doc: free>] 271 patients: [] doctors: [<Doc: time left=14>, <Doc: free>, <Doc: free>] 272 patients: [] doctors: [<Doc: time left=13>, <Doc: free>, <Doc: free>] 273 patients: [] doctors: [<Doc: time left=12>, <Doc: free>, <Doc: free>] 274 patients: [] doctors: [<Doc: time left=11>, <Doc: free>, <Doc: free>] 275 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 276 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 277 patients: [] doctors: [<Doc: time left=8>, <Doc: free>, <Doc: free>] 278 patients: [] doctors: [<Doc: time left=7>, <Doc: free>, <Doc: free>] 279 patients: [] doctors: [<Doc: time left=6>, <Doc: free>, <Doc: free>] 280 patients: [] doctors: [<Doc: time left=5>, <Doc: free>, <Doc: free>] 281 patients: [] doctors: [<Doc: time left=4>, <Doc: free>, <Doc: free>] 282 patients: [] doctors: [<Doc: time left=3>, <Doc: free>, <Doc: free>] 283 patient arrived! patients: [] doctors: [<Doc: time left=2>, <Doc: time left=22>, <Doc: free>] 284 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=21>, <Doc: free>] 285 patients: [] doctors: [<Doc: free>, <Doc: time left=20>, <Doc: free>] 286 patients: [] doctors: [<Doc: free>, <Doc: time left=19>, <Doc: free>] 287 patients: [] doctors: [<Doc: free>, <Doc: time left=18>, <Doc: free>] 288 patients: [] doctors: [<Doc: free>, <Doc: time left=17>, <Doc: free>] 289 patients: [] doctors: [<Doc: free>, <Doc: time left=16>, <Doc: free>] 290 patients: [] doctors: [<Doc: free>, <Doc: time left=15>, <Doc: free>] 291 patients: [] doctors: [<Doc: free>, <Doc: time left=14>, <Doc: free>] 292 patients: [] doctors: [<Doc: free>, <Doc: time left=13>, <Doc: free>] 293 patients: [] doctors: [<Doc: free>, <Doc: time left=12>, <Doc: free>] 294 patients: [] doctors: [<Doc: free>, <Doc: time left=11>, <Doc: free>] 295 patients: [] doctors: [<Doc: free>, <Doc: time left=10>, <Doc: free>] 296 patients: [] doctors: [<Doc: free>, <Doc: time left=9>, <Doc: free>] 297 patients: [] doctors: [<Doc: free>, <Doc: time left=8>, <Doc: free>] 298 patients: [] doctors: [<Doc: free>, <Doc: time left=7>, <Doc: free>] 299 patients: [] doctors: [<Doc: free>, <Doc: time left=6>, <Doc: free>] 300 patient arrived! patients: [] doctors: [<Doc: time left=18>, <Doc: time left=5>, <Doc: free>] 301 patients: [] doctors: [<Doc: time left=17>, <Doc: time left=4>, <Doc: free>] 302 patients: [] doctors: [<Doc: time left=16>, <Doc: time left=3>, <Doc: free>] 303 patients: [] doctors: [<Doc: time left=15>, <Doc: time left=2>, <Doc: free>] 304 patients: [] doctors: [<Doc: time left=14>, <Doc: time left=1>, <Doc: free>] 305 patients: [] doctors: [<Doc: time left=13>, <Doc: free>, <Doc: free>] 306 patients: [] doctors: [<Doc: time left=12>, <Doc: free>, <Doc: free>] 307 patients: [] doctors: [<Doc: time left=11>, <Doc: free>, <Doc: free>] 308 patients: [] doctors: [<Doc: time left=10>, <Doc: free>, <Doc: free>] 309 patients: [] doctors: [<Doc: time left=9>, <Doc: free>, <Doc: free>] 310 patients: [] doctors: [<Doc: time left=8>, <Doc: free>, <Doc: free>] 311 patients: [] doctors: [<Doc: time left=7>, <Doc: free>, <Doc: free>] 312 patients: [] doctors: [<Doc: time left=6>, <Doc: free>, <Doc: free>] 313 patients: [] doctors: [<Doc: time left=5>, <Doc: free>, <Doc: free>] 314 patients: [] doctors: [<Doc: time left=4>, <Doc: free>, <Doc: free>] 315 patients: [] doctors: [<Doc: time left=3>, <Doc: free>, <Doc: free>] 316 patients: [] doctors: [<Doc: time left=2>, <Doc: free>, <Doc: free>] 317 patients: [] doctors: [<Doc: time left=1>, <Doc: free>, <Doc: free>] 318 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 319 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 320 patient arrived! patients: [] doctors: [<Doc: time left=22>, <Doc: free>, <Doc: free>] 321 patients: [] doctors: [<Doc: time left=21>, <Doc: free>, <Doc: free>] 322 patients: [] doctors: [<Doc: time left=20>, <Doc: free>, <Doc: free>] 323 patients: [] doctors: [<Doc: time left=19>, <Doc: free>, <Doc: free>] 324 patient arrived! patients: [] doctors: [<Doc: time left=18>, <Doc: time left=14>, <Doc: free>] 325 patients: [] doctors: [<Doc: time left=17>, <Doc: time left=13>, <Doc: free>] 326 patients: [] doctors: [<Doc: time left=16>, <Doc: time left=12>, <Doc: free>] 327 patients: [] doctors: [<Doc: time left=15>, <Doc: time left=11>, <Doc: free>] 328 patients: [] doctors: [<Doc: time left=14>, <Doc: time left=10>, <Doc: free>] 329 patients: [] doctors: [<Doc: time left=13>, <Doc: time left=9>, <Doc: free>] 330 patients: [] doctors: [<Doc: time left=12>, <Doc: time left=8>, <Doc: free>] 331 patients: [] doctors: [<Doc: time left=11>, <Doc: time left=7>, <Doc: free>] 332 patients: [] doctors: [<Doc: time left=10>, <Doc: time left=6>, <Doc: free>] 333 patients: [] doctors: [<Doc: time left=9>, <Doc: time left=5>, <Doc: free>] 334 patients: [] doctors: [<Doc: time left=8>, <Doc: time left=4>, <Doc: free>] 335 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=3>, <Doc: free>] 336 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=2>, <Doc: free>] 337 patient arrived! patients: [] doctors: [<Doc: time left=5>, <Doc: time left=1>, <Doc: time left=24>] 338 patients: [] doctors: [<Doc: time left=4>, <Doc: free>, <Doc: time left=23>] 339 patients: [] doctors: [<Doc: time left=3>, <Doc: free>, <Doc: time left=22>] 340 patient arrived! patients: [] doctors: [<Doc: time left=2>, <Doc: time left=21>, <Doc: time left=21>] 341 patient arrived! patients: [<patient: wait=1>] doctors: [<Doc: time left=1>, <Doc: time left=20>, <Doc: time left=20>] 342 patients: [] doctors: [<Doc: time left=25>, <Doc: time left=19>, <Doc: time left=19>] 343 patients: [] doctors: [<Doc: time left=24>, <Doc: time left=18>, <Doc: time left=18>] 344 patients: [] doctors: [<Doc: time left=23>, <Doc: time left=17>, <Doc: time left=17>] 345 patients: [] doctors: [<Doc: time left=22>, <Doc: time left=16>, <Doc: time left=16>] 346 patient arrived! patients: [<patient: wait=1>] doctors: [<Doc: time left=21>, <Doc: time left=15>, <Doc: time left=15>] 347 patients: [<patient: wait=2>] doctors: [<Doc: time left=20>, <Doc: time left=14>, <Doc: time left=14>] 348 patients: [<patient: wait=3>] doctors: [<Doc: time left=19>, <Doc: time left=13>, <Doc: time left=13>] 349 patients: [<patient: wait=4>] doctors: [<Doc: time left=18>, <Doc: time left=12>, <Doc: time left=12>] 350 patients: [<patient: wait=5>] doctors: [<Doc: time left=17>, <Doc: time left=11>, <Doc: time left=11>] 351 patients: [<patient: wait=6>] doctors: [<Doc: time left=16>, <Doc: time left=10>, <Doc: time left=10>] 352 patients: [<patient: wait=7>] doctors: [<Doc: time left=15>, <Doc: time left=9>, <Doc: time left=9>] 353 patients: [<patient: wait=8>] doctors: [<Doc: time left=14>, <Doc: time left=8>, <Doc: time left=8>] 354 patients: [<patient: wait=9>] doctors: [<Doc: time left=13>, <Doc: time left=7>, <Doc: time left=7>] 355 patients: [<patient: wait=10>] doctors: [<Doc: time left=12>, <Doc: time left=6>, <Doc: time left=6>] 356 patients: [<patient: wait=11>] doctors: [<Doc: time left=11>, <Doc: time left=5>, <Doc: time left=5>] 357 patients: [<patient: wait=12>] doctors: [<Doc: time left=10>, <Doc: time left=4>, <Doc: time left=4>] 358 patients: [<patient: wait=13>] doctors: [<Doc: time left=9>, <Doc: time left=3>, <Doc: time left=3>] 359 patients: [<patient: wait=14>] doctors: [<Doc: time left=8>, <Doc: time left=2>, <Doc: time left=2>] 360 patients: [<patient: wait=15>] doctors: [<Doc: time left=7>, <Doc: time left=1>, <Doc: time left=1>] 361 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=20>, <Doc: free>] 362 patients: [] doctors: [<Doc: time left=5>, <Doc: time left=19>, <Doc: free>] 363 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=18>, <Doc: free>] 364 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=17>, <Doc: free>] 365 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=16>, <Doc: free>] 366 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=15>, <Doc: free>] 367 patients: [] doctors: [<Doc: free>, <Doc: time left=14>, <Doc: free>] 368 patient arrived! patients: [] doctors: [<Doc: time left=12>, <Doc: time left=13>, <Doc: free>] 369 patient arrived! patients: [] doctors: [<Doc: time left=11>, <Doc: time left=12>, <Doc: time left=16>] 370 patients: [] doctors: [<Doc: time left=10>, <Doc: time left=11>, <Doc: time left=15>] 371 patients: [] doctors: [<Doc: time left=9>, <Doc: time left=10>, <Doc: time left=14>] 372 patients: [] doctors: [<Doc: time left=8>, <Doc: time left=9>, <Doc: time left=13>] 373 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=8>, <Doc: time left=12>] 374 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=7>, <Doc: time left=11>] 375 patients: [] doctors: [<Doc: time left=5>, <Doc: time left=6>, <Doc: time left=10>] 376 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=5>, <Doc: time left=9>] 377 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=4>, <Doc: time left=8>] 378 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=3>, <Doc: time left=7>] 379 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=2>, <Doc: time left=6>] 380 patients: [] doctors: [<Doc: free>, <Doc: time left=1>, <Doc: time left=5>] 381 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=4>] 382 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=3>] 383 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=2>] 384 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: time left=1>] 385 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 386 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 387 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 388 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 389 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 390 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 391 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 392 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 393 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 394 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 395 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 396 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 397 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 398 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 399 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 400 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 401 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 402 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 403 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 404 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 405 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 406 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 407 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 408 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>] 409 patient arrived! patients: [] doctors: [<Doc: time left=24>, <Doc: free>, <Doc: free>] 410 patients: [] doctors: [<Doc: time left=23>, <Doc: free>, <Doc: free>] 411 patients: [] doctors: [<Doc: time left=22>, <Doc: free>, <Doc: free>] 412 patients: [] doctors: [<Doc: time left=21>, <Doc: free>, <Doc: free>] 413 patients: [] doctors: [<Doc: time left=20>, <Doc: free>, <Doc: free>] 414 patients: [] doctors: [<Doc: time left=19>, <Doc: free>, <Doc: free>] 415 patients: [] doctors: [<Doc: time left=18>, <Doc: free>, <Doc: free>] 416 patients: [] doctors: [<Doc: time left=17>, <Doc: free>, <Doc: free>] 417 patients: [] doctors: [<Doc: time left=16>, <Doc: free>, <Doc: free>] 418 patients: [] doctors: [<Doc: time left=15>, <Doc: free>, <Doc: free>] 419 patient arrived! patients: [] doctors: [<Doc: time left=14>, <Doc: time left=25>, <Doc: free>] 420 patients: [] doctors: [<Doc: time left=13>, <Doc: time left=24>, <Doc: free>] 421 patients: [] doctors: [<Doc: time left=12>, <Doc: time left=23>, <Doc: free>] 422 patients: [] doctors: [<Doc: time left=11>, <Doc: time left=22>, <Doc: free>] 423 patients: [] doctors: [<Doc: time left=10>, <Doc: time left=21>, <Doc: free>] 424 patients: [] doctors: [<Doc: time left=9>, <Doc: time left=20>, <Doc: free>] 425 patients: [] doctors: [<Doc: time left=8>, <Doc: time left=19>, <Doc: free>] 426 patients: [] doctors: [<Doc: time left=7>, <Doc: time left=18>, <Doc: free>] 427 patients: [] doctors: [<Doc: time left=6>, <Doc: time left=17>, <Doc: free>] 428 patients: [] doctors: [<Doc: time left=5>, <Doc: time left=16>, <Doc: free>] 429 patients: [] doctors: [<Doc: time left=4>, <Doc: time left=15>, <Doc: free>] 430 patients: [] doctors: [<Doc: time left=3>, <Doc: time left=14>, <Doc: free>] 431 patients: [] doctors: [<Doc: time left=2>, <Doc: time left=13>, <Doc: free>] 432 patients: [] doctors: [<Doc: time left=1>, <Doc: time left=12>, <Doc: free>] 433 patients: [] doctors: [<Doc: free>, <Doc: time left=11>, <Doc: free>] 434 patients: [] doctors: [<Doc: free>, <Doc: time left=10>, <Doc: free>] 435 patients: [] doctors: [<Doc: free>, <Doc: time left=9>, <Doc: free>] 436 patients: [] doctors: [<Doc: free>, <Doc: time left=8>, <Doc: free>] 437 patients: [] doctors: [<Doc: free>, <Doc: time left=7>, <Doc: free>] 438 patients: [] doctors: [<Doc: free>, <Doc: time left=6>, <Doc: free>] 439 patients: [] doctors: [<Doc: free>, <Doc: time left=5>, <Doc: free>] 440 patients: [] doctors: [<Doc: free>, <Doc: time left=4>, <Doc: free>] 441 patients: [] doctors: [<Doc: free>, <Doc: time left=3>, <Doc: free>] 442 patients: [] doctors: [<Doc: free>, <Doc: time left=2>, <Doc: free>] 443 patients: [] doctors: [<Doc: free>, <Doc: time left=1>, <Doc: free>] 444 patients: [] doctors: [<Doc: free>, <Doc: free>, <Doc: free>]
{'avg_wait': 5.666666666666667, 'longest_wait': 15, 'mins_past_close': 25, 'n_patients': 30, 'n_waiters': 3, 'wait_times': [0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 15, 0, 0, 0, 0]}
data = [run_simulation() for _ in range(2000)]
patient_counts = [x['n_patients'] for x in data]
avg_waits = [x.get('avg_wait', 0) for x in data]
longest_waits = [x['longest_wait'] for x in data]
waiter_counts = [x['n_waiters'] for x in data]
overtimes = [x['mins_past_close'] for x in data if x['mins_past_close'] > 0]
def get_percentiles(data):
return np.percentile(data, [25, 50, 75])
for stat in [patient_counts, avg_waits, longest_waits,
waiter_counts, overtimes]:
print(get_percentiles(stat))
[ 36. 40. 44.] [ 3.5 5. 6.5] [ 6. 9. 13.] [ 3. 6. 10.] [ 7. 12. 17.25]
plt.rcParams['figure.figsize'] = 7, 9
plt.subplot(3, 2, 1)
plt.hist(patient_counts, color='black')
plt.xlabel('number of patients')
plt.title('Patient Counts')
plt.subplot(3, 2, 2)
plt.hist(avg_waits, color='black')
plt.xlabel('minutes')
plt.title('Average Waits')
plt.subplot(3, 2, 3)
plt.hist(longest_waits, color='black')
plt.xlabel('minutes')
plt.title('Longest Waits')
plt.subplot(3, 2, 4)
plt.hist(waiter_counts, color='black')
plt.xlabel('number of waiters')
plt.title('Number of people who waited')
plt.subplot(3, 2, 5)
plt.hist(overtimes, color='black')
plt.xlabel('minutes')
plt.title('Number of minutes Overtime')
plt.tight_layout();
# plt.savefig('stats_2.png', format='png', dpi=400)
2. Consider two coins, C1 and C2, with the following characteristics: Pr(heads|C1) = 0.4 and Pr(heads|C2) = 0.6. Choose one of the coins at random and imagine spinning it repeatedly.
1. Given that the first two spins from the chosen coin are tails, what is the expected number of additional spins until a heads shows up?
2. How about if the first three spins are all tails?
x = np.linspace(0, 1, 100)
prior = beta.pdf(x, 3, 3)
likelihood = binom.cdf(3, 10, x)
posterior = likelihood * prior
lw = 4
color = 'black'
plt.rcParams['figure.figsize'] = 8, 6
plt.style.use('ggplot')
plt.subplot(2, 1, 1)
plt.plot(x, prior, linewidth=lw, color=color)
plt.title('Prior')
plt.xlabel(r'$\theta$', labelpad=.01)
plt.ylabel(r'$p(\theta)$')
plt.yticks([])
plt.grid(b=False)
plt.subplot(2, 1, 2)
plt.plot(x, posterior, linewidth=lw, color=color)
plt.title('Posterior')
plt.xlabel(r'$\theta$')
plt.ylabel(r'$p(\theta|y)$')
plt.yticks([])
plt.grid(b=False)
plt.tight_layout();
plt.savefig('coin_flip_densities.png', format='png', dpi=400)
df = pd.read_csv('data/accident_data.csv', index_col='year')
df
accidents | exposures | |
---|---|---|
year | ||
1976 | 24 | 3864 |
1977 | 25 | 4300 |
1978 | 31 | 5026 |
1979 | 31 | 5482 |
1980 | 22 | 5815 |
1981 | 21 | 6034 |
1982 | 26 | 5876 |
1983 | 20 | 6224 |
1984 | 16 | 7434 |
1985 | 22 | 7106 |
# check if the accident rate is increasing or decreasing
rates = df['accidents'] / df['exposures']
x = df.index.tolist()
# plot it
plt.scatter(x, rates, s=70, color='black')
plt.xlim([1975, 1986])
plt.ylim([0, 0.01])
plt.xlabel('Year')
plt.ylabel('Accident Rate')
plt.title('Accident Rates by Year')
plt.show();
The following table gives the number of fatal accidents on scheduled airline flights, and exposure (in 100s of millions of passenger miles), per year, over a ten-year period.
Assume that the numbers of fatal accidents in each year are independent with a Poisson(θ) distribution. Assign a $Gamma(\alpha,\beta)$ prior for $\theta$, selecting reasonable values for $\alpha$ and $\beta$. Using simulation, compute a 95% predictive interval for the number of fatal acci- dents in 1986.
For t = 1,...,10, let yt denote the number of fatal accidents in year 1975 + t. Assume the yt are independent $Poisson(xt\theta)$, where xt denotes the exposure (in 100s of millions of passenger miles) in year 1975 + t. Assign a $Gamma(\alpha,\beta)$ prior for $\theta$; you will probably want different hyperparameter values than in part (1). Give a 95% predictive interval for the number of fatal accidents in 1986 under the assumption that 7.6 × 1011 passenger miles are flown that year.
$\theta$ ~ $Gamma(\alpha, \beta)$
$y_{i}$ ~ $Poisson(\theta)$
class BayesianPoissonModel:
'''Spits out data simulated values from a bayesian model'''
def __init__(self, initial_alpha, initial_beta, n_exposures=1):
self.alpha = initial_alpha
self.beta = initial_beta
self.n_exposures = n_exposures
def update(self, data):
'''
updates p(theta) to be p(theta|data)
'''
for y, n in data:
self.alpha += y
self.beta += n
def simulate(self, n):
theta = np.random.gamma(self.alpha, self.beta**-1, n)
return np.random.poisson(theta * self.n_exposures)
def __repr__(self):
return '<model: alpha={:.2f}, beta={:.2f}>'.format(self.alpha, self.beta)
def __str__(self):
return self.__repr__()
We'll disregard number of exposures per year and only model the number of accidents
$\beta$ will be 1, as in a period of 1 year
accidents_data = np.c_[df['accidents'].values, np.ones(df.shape[0])]
model = BayesianPoissonModel(df['accidents'].mean(), 1)
model.update(accidents_data)
simulated = model.simulate(n=100000)
print(model)
print('95% predictive interval: [{}, {}]'.format(*np.percentile(simulated, (2.5, 97.5)).astype(int)))
# plot it
plt.rcParams['figure.figsize'] = 7, 4
plt.hist(simulated, bins=20, color='black')
plt.xlabel('Number of Accidents')
plt.ylabel('Count')
plt.title('Simulated Number of Fatal Accidents')
plt.show();
<model: alpha=261.80, beta=11.00> 95% predictive interval: [14, 34]
Incorporate number of exposures into the models to "weight" each accident count appropriately
n_exposures = 7600
# use the means as a prior
a = df['accidents'].mean()
b = df['exposures'].mean()
model_full = BayesianPoissonModel(a, b, n_exposures=n_exposures)
model_full.update(df[['accidents', 'exposures']].values)
simulated = model_full.simulate(n=100000)
print(model_full)
print('95% predictive interval: [{}, {}]'.format(*np.percentile(simulated, (2.5, 97.5)).astype(int)))
# plot it
plt.rcParams['figure.figsize'] = 7, 4
plt.hist(simulated, bins=20, color='black')
plt.xlabel('Number of Accidents')
plt.ylabel('count')
plt.title('Simulated Number of Fatal Accidents')
plt.show;
<model: alpha=261.80, beta=62877.10> 95% predictive interval: [21, 44]
An experiment was performed to estimate the effect of beta blockers on mortality of cardiac patients. A group of patients were randomly assigned to treatment and control groups: out of 282 patients receiving the control, 37 died; out of 278 receiving the treatment, 25 died. Assume that the outcomes are independent and binomially distributed, with probabilities of death $p_0$ and $p_1$ under the control and treatment, respectively.
Set up a noninformative prior distribution on ($p_0$, $p_1$) and obtain posterior simulations. Summarize the posterior distribution for the odds ratio,
# draw 1000 p_0's and p_1's from their respective posteriors
n_samples = 1000
p_0 = np.random.beta(37 + 1, 282 - 37 + 1, n_samples) # control
p_1 = np.random.beta(25 + 1, 278 - 25 + 1, n_samples) # treatment
# compute each simulated odds ratio
odds_samples = (p_1 / (1 - p_1)) / (p_0 / (1 - p_0))
lower, upper = np.percentile(odds_samples, [2.5, 97.5])
print('95% predictive interval: [{:.3f}, {:.3f}]'.format(lower, upper))
# plot it
plt.hist(odds_samples, bins=15, color='black')
plt.xlabel('Odds')
plt.ylabel('frequency')
plt.title('Simulated odds ratios of $p_1$ and $p_0$', y=1.03)
# plt.savefig('p_7_plot.png', format='png', dpi=400)
plt.show()
95% predictive interval: [0.403, 1.113]