Monte Carlo simulation of the Monty Hall Problem in Ruby and Python

Reading Jeff Atwood’s post The Problem of the Unfinished Game, reminded me of a similar problem. The Monty Hall Problem is a well known probability puzzle that has tricked many people. In fact, if you are not familiar with it already, chances are that you’ll get it wrong. And you would be in good company along with many mathematicians and physicists, including the great mathematician, Paul Erdos. This puzzle is loosely based on the television show Let’s Make a Deal, and is equivalent to some much older puzzles you may be familiar with (e.g. the three prisoners problem). In its simplest form, it asks the following question:

Suppose you’re on a game show, and you’re given the choice of three doors: Behind one door is a car; behind the others, goats. You pick a door, say No. 1, and the host, who knows what’s behind the doors, opens another door, say No. 3, which has a goat. He then says to you, “Do you want to pick door No. 2?” Is it to your advantage to switch your choice?

This definition of the problem is admittedly ambiguous. Thankfully Wikipedia points us towards a more exact definition:

Suppose you’re on a game show and you’re given the choice of three doors. Behind one door is a car; behind the others, goats [that is, booby prizes]. The car and the goats were placed randomly behind the doors before the show. The rules of the game show are as follows: After you have chosen a door, the door remains closed for the time being. The game show host, Monty Hall, who knows what is behind the doors, now has to open one of the two remaining doors, and the door he opens must have a goat behind it. If both remaining doors have goats behind them, he chooses one randomly. After Monty Hall opens a door with a goat, he will ask you to decide whether you want to stay with your first choice or to switch to the last remaining door. Imagine that you chose Door 1 and the host opens Door 3, which has a goat. He then asks you “Do you want to switch to Door Number 2?” Is it to your advantage to change your choice?

Think about it for a moment, then read on. To answer this question, most people will try to determine which of the two possible outcomes has a higher probability. Problems arise when trying to correctly calculate the probability of these two events though. There are two closed doors and the car could be behind either of them. Hence, most people’s “common sense” and psychology leads them to believe that there is a 50% chance that the car is behind the initially selected door, and 50% that it’s behind the other closed door that was offered up by Monty. Initially it would seem that switching or staying with the first choice doesn’t really make a difference.

Unfortunately that’s not the right answer. The correct answer is that there is a two out of three chance of winning by switching to the other door; so switching is always to your advantage. This result is considered to be a paradox because it’s very counterintuitive to the way that many people think. It is in fact so counterintuitive that most people will argue with you in an attempt to convince you otherwise. I invite you to check out the Wikipedia entry on the problem/paradox, to read a step-by-step explanation with figures about why switching gives you about 66.7% chance of winning the car and why staying with the initial choice gives you only a 33.3% success rate.

When you make your first choice your probability of winning the car is only 1/3. If you decide to switch, you will win only if the first choice you made was wrong. And since your first choice came with a 2 out of 3 chance of picking a goat, switching will then (logically) give you 2/3 chance of winning. Another easy way to come to intuitively accept this surprising result, is to wildly exaggerate the terms of the problem. If there were a billion doors, you picked one, and then Monty proceeded to open up all the remaining doors but one, we’d have a situation where it would be extremely unlikely that you picked the right door at the beginning, while it would be extremely likely that the remaining door was the one that was concealing the car.

Even after reading several explanations and aids to understand these results, there are still people who are skeptical or refuse to believe them. Let’s verify the outcome with a simulation.

What you find below is a quick Ruby script that I wrote to run a Monte Carlo Simulation of the Monty Hall problem/paradox. It runs the game a million times and then measures how many times the player won by sticking with their first choice, and how many times switching would have led to winning the car.

#!/usr/bin/env ruby -w

# Monte Carlo simulation for the Monty Hall Problem:
# http://en.wikipedia.org/wiki/Monty_Hall_problem

=begin
When using a Ruby version older than 1.8.7
define the following two methods:

  class Array
    def shuffle
      self.sort_by { rand }
    end
    
    def choice
      self.shuffle.first
    end
  end
=end

# Utility class for the simulation of a single Monty Hall game.
class MontyHall
  def initialize
    @doors = ['car', 'goat', 'goat'].shuffle
  end

  # Return a number representing the player's first choice.
  def pick_door
    return rand(3)
  end

  # Return the index of the door opened by the host.
  # This cannot represent a door hiding a car or the player's chosen door.
  def reveal_door(pick)
    available_doors = [0, 1, 2]
    available_doors.delete(pick)
    available_doors.delete(@doors.index('car'))
    return available_doors.choice
  end

  # Return true if the player won by staying
  # with their first choice, false otherwise.
  def staying_wins?(pick)
    won?(pick)
  end

  # Return true if the player won by switching, false otherwise.
  def switching_wins?(pick, open_door)
    switched_pick = ([0, 1, 2] - [open_door, pick]).first
    won?(switched_pick)
  end

  private

  # Return true if the player's final pick hides a car, false otherwise.
  def won?(pick)
    @doors[pick] == 'car'
  end
end

if __FILE__ == $0
  ITERATIONS = (ARGV.shift || 1_000_000).to_i
  staying = 0
  switching = 0

  ITERATIONS.times do
    mh = MontyHall.new
    picked = mh.pick_door
    revealed = mh.reveal_door(picked)
    staying += 1 if mh.staying_wins?(picked)
    switching += 1 if mh.switching_wins?(picked, revealed)
  end

  staying_rate = (staying.to_f / ITERATIONS) * 100
  switching_rate = (switching.to_f / ITERATIONS) * 100

  puts "Staying: #{staying_rate}%."
  puts "Switching: #{switching_rate}%."
end

And here is an “equivalent” version I wrote in Python:

#!/usr/bin/env python
"""
Monte Carlo simulation for the Monty Hall Problem:
http://en.wikipedia.org/wiki/Monty_Hall_problem.
"""
import sys
from random import randrange, shuffle, choice

DOORS = ['car', 'goat', 'goat']

def pick_door():
    """Return a number representing the player's first choice."""
    return randrange(3)

def reveal_door(pick):
    """Return the index of the door opened by the host.
    This cannot be a door hiding a car or the player's chosen door.
    """
    all_doors = set([0, 1, 2])
    unavailable_doors = set([DOORS.index('car'), pick])
    available_doors = list(all_doors - unavailable_doors)
    return choice(available_doors)

def staying_wins(pick):
    """Return True if the player won by staying
    with their first choice, False otherwise.
    """
    return won(pick)

def switching_wins(pick, open_door):
    """Return True if the player won by switching,
    False otherwise.
    """
    other_doors = set([pick, open_door])
    switched_pick = (set([0, 1, 2]) - other_doors).pop()
    return won(switched_pick)

def won(pick):
    """Return True if the player's final pick hides a car,
    False otherwise.
    """
    return (DOORS[pick] == 'car')

def main(iterations=1000000):
    """Run the main simulation as many
    times as specified by the function argument.
    """
    shuffle(DOORS)

    switching = 0
    staying = 0

    for dummy in xrange(iterations):
        picked = pick_door()
        revealed = reveal_door(picked)
        if staying_wins(picked):
            staying += 1
        if switching_wins(picked, revealed):
            switching += 1

    staying_rate = (float(staying) / iterations) * 100
    switching_rate = (float(switching) / iterations) * 100

    print "Staying: %f%%" % staying_rate
    print "Switching: %f%%" % switching_rate

if __name__ == "__main__":
    if len(sys.argv) == 2:
        main(int(sys.argv[1]))
    else:
        main()

Even if you are not familiar with Ruby or Python, you may be able to understand what’s going on here. The main body of the program emulates the game and keeps track of the number of victories when the player sticks with their initial choice, and when they switch. Notice that this code intentionally tries not to be clever, in order not to annoy “skeptical” people.

There are many points in the code where correct assumptions about the problem would lead us to code that is faster and much more compact. For example, if the player wins a given game by sticking with his first answer, it’s obvious that switching would have made him lose. We could just calculate the difference between 100 and the success rate of staying with the first choice, and we’d obtain the success rate for switching. But here we are trying to simulate the problem as faithfully as possible and abstract as little as necessary.

As always with Monte Carlo Simulations, the outcome is slightly variable during each run since it depends on random input; but by the law of large numbers, it will very slowly converge to the expected values (despite the pseudo-randomness used here). For example, when I executed the code above for the first time on my machine, I obtained the following:

Staying: 33.382%.
Switching: 66.618%.

The results of this simulation should be enough to convince you that the theoretical results are actually true; we are easily fooled, and the mathematicians who got it right were not making stuff up.

Happy New Year to my readers, I wish you all the best for a happy, successful 2009!

How to parse decimal numbers within a string

INPUT: a string containing decimal numbers.

OUTPUT: an array containing all the decimal numbers within the given string.

You can accomplish this task very quickly with the String#scan method and the right Regular Expression (regex).

Given a string s, you can use:

numbers = s.scan /[-+]?\d*\.?\d+/

numbers will be an array whose elements are the decimal numbers within the string s. Note how the regex considers the possible + or – signs in front of the numbers.

If you also wish to match floating point numbers with exponents (scientific notation, e.g. 2.54.e-07), then use the following:

numbers = s.scan /[-+]?\d*\.?\d+([eE][-+]?\d+)?/ 

Solving mathematical problems with Ruby

A couple of weeks ago I’ve discovered a very interesting and addictive website MathsChallenge.net . They have a section called Project Euler where participants can try to resolve more than 100 mathematical problems. The scoring system is really well done and the discussion forum for each problem is accessible only after having solved the given problem.

It’s interesting to see how other people try to figure out a solution and what languages they employ to solve the puzzling problems. I even discovered rare languages like J or K through this site. It’s like a martial arts contest where participants can use their favorite martial art, only this is for hackers.

I’ve only spent a few hours over a couple of days on the site, enough time to solve 26 problems and bring me somewhere around 370th place out of more than 1300 active participants in the contest. I didn’t spend time solving problems during the last ten days, and in fact I am now 418th. I am going to work on more problems when I’ll have some time. I love this kind of stuff, and it’s not just for the enjoyable competition. It doesn’t really matter if you are first or last, especially if you know how to solve most of the problems but you are short of time like me. What’s important is to have fun and benefit from the learning process.

The programs must solve the given problems within 60 seconds (think of it like Gone in 60 Seconds), which means that most brute force methods won’t work.

This is particularly true for Ruby. In fact I have seen some people with C or Assembler resolve a problem with brute force algorithms, where Ruby would take way longer than a minute. There is no harm in this though, it’s a good exercise in controlling the algorithms’ complexity, facing problems that you don’t encounter too often in real world programming.

I use Ruby 1.8.4 so I even had to rewrite my own mathn library because prime computation and gcd are exceptionally slow (1.9 solves the problem though).

I invite you to join me and bring up Ruby statistics on MathsChallenge. And above all, have lots of fun…