Proving the Monty Hall Problem, by simulating 1 million runs in python



New to the monty hall problem?

In short, the problem is:

In the game show, there are 3 doors. 1 has a car, 2 have sheep. You get whatever is in the door you picked.

The host knows which door has the car, and which doors have sheep.

You choose a door. The host opens another door and shows you that it has sheep.

Now, you are given a choice:

  1. Stick to the door you picked
  2. Switch to the other remaining door

What has a higher probability of winning the car, sticking or changing?

Need more details? read it its wikipedia article


The math says that before the reveal, your probability of winning was 33%. But after the reveal, if you switch your door - your probability of winning becomes 66%.

This is absolutely mind bending.

I looked at the math. But was not really convinced.

You know what's better than modelling the problem in math, solving it, and coming up with the final probabilities? Write the simulation in code, run it a million times, and directly look at the resultant probabilities!

So I programmed a simulation, ran it 1 million times, and looked at the results. Now this is a method that I automatically trust.

I can be bad at math, but I cannot be bad at directly looking at the results with my own eyes.

Result: it is true. Changing the door gives a higher chance of winning.

normal_probability:  0.334064
stick_door_probability:  0.334064
change_door_probability:  0.665936

The entire program, if you want it:

# Usual warning: code was written in 2-3 minutes.
# It is slow (no numpy)
# It is not clean (again, wrote it in 2-3 minutes)
# If you write better code, email it to me or write your own blog post. That's it. Enjoy.

import random
iterations = 1_000_000
# generate door simulations
simulations = np.full((iterations, 3), False)
# generate column of true
index_of_true = np.random.randint(3, size=iterations)
# set true values
simulations[range(iterations), index_of_true] = True

simulations_array = simulations.tolist()

normal_results = []
change_door_results = []
stick_door_results = []
for simulation in simulations_array:
    # pick a random door
    # check if door is correct or wrong
    chosen_door = random.randint(0, 2)

    # flippable door situration
    # eliminate a false door other than chosen door
    #   pick a closed door

    all_doors = [0,1,2]
    non_chosen_doors = all_doors

    sheep_doors_in_non_chosen_doors = []
    for door_index in non_chosen_doors:
        if not simulation[door_index]:

    # choose a random sheep door to eliminate
    door_eliminated = random.choice(sheep_doors_in_non_chosen_doors)

    # don't change door, and check if guess was right

    # change door, and check if guess was right
    all_doors = [0, 1, 2]
    other_door_index = all_doors[0]

    changed_chosen_door = other_door_index

# calculate probability of winning
normal_probability = sum(normal_results) / len(normal_results)
stick_door_probability = sum(stick_door_results) / len(stick_door_results)
change_door_probability = sum(change_door_results) / len(change_door_results)

print(f"normal_probability: {normal_probability}"")
print(f"stick_door_probability: {stick_door_probability_}")
print(f"change_door_probability: {change_door_probability}")

That's the fun of programming. You can quickly simulate most probability based problems and see things for yourself.