Tuesday, November 29, 2011

A Birthday Simulator

Earlier today I was reading about the Tuesday Birthday Problem (which curiously doesn't seem to have its own entry on Wikipedia.. maybe it is known under a different name?) and although I was convinced by the argument, I thought that a little simulation would help deepen my understanding of this strange paradox (or at least make it a little more intuitive). The problem I had is how to represent, in a clear way, some a priori knowledge (namely, the fact that one of the children is a son born on a Tuesday) in a numerical simulation?

Since directly modeling the conditional distribution wouldn't be trivial, an easier way to do it is by using rejection sampling: iterate over a set of randomly generated family configurations, and reject those that do not match the given fact, i.e. those not containing at least a son born on a Tuesday. From the configurations that passed the test, the proportion of those having the other child also a son (born on whatever day), should yield the answer (which of course is not 1/2, as intuition first strongly suggests):

from __future__ import division
from random import *

n = 0
n_times_other_child_is_son = 0
while n < 100000:
    child1 = choice('MF') + str(randint(0, 6))
    child2 = choice('MF') + str(randint(0, 6))
    children = child1 + child2
    if 'M2' not in children: continue
    if children[0::2] == 'MM':
        n_times_other_child_is_son += 1
    n += 1

print n_times_other_child_is_son / n # should be close to 13/27