Problem

Suppose you’re on a game show, and you’re given the choice of three doors,


Behind one door is a car; behind the two others, goats. You pick a door, say No. 1, and the host of the show opens another door, say No. 3, which has a goat.


He then says to you, “Do you want to pick door No. 2?”.

Question: What would you do?

Is it to your advantage to switch your choice from door 1 to door 2? Is it to your advantage, in the long run, for a large number of game-tries, to switch to the other door?

Now whatever your answer is, I want you to check/prove your answer by a Monte Carlo simulation of this problem. Make a plot of your simulation for $\rm{ngames}=100000$ repeat of this game, that shows, in the long run, on average, what the probability of winning this game is if you switch your choice, and what is the probability of winning, if you do not switch to the other door.

Solution

As you see in the figure, although you may initially win by not switching your choice, in the long run, on average, you will lose, if you don’t switch your choice.

MATLAB

Here is an example MATLAB implementation, monteGame.m,

%rng('shuffle');
rng(1235);
nExperiments=100000; % the number of times the experiment is simulated
TotalWinsWithoutSwitch = 0;
TotalWinsWithSwitch = 0;
AllDoors = [1,2,3];
averageNumberOfWinsWithoutSwitch = zeros(nExperiments,1);
averageNumberOfWinsWithSwitch = zeros(nExperiments,1);

for iExperiment = 1:nExperiments
        
    if mod(iExperiment,10000)==0; disp(['simulating experiment number ', num2str(iExperiment)]); end
    
    ContainsCar = zeros(3,1); % first create the three doors, assuming none contains the car, hence all are zero valued.
    doorWithCar = randi([1,3],1,1); % Now put the car behind one of the three doors randomly.
    myChoice = randi([1,3],1,1); % Now make a choice of door randomly.
    
    excludedDoor = AllDoors(AllDoors~=doorWithCar & AllDoors~=myChoice); % find which door should be opened by the host
    if length(excludedDoor)>1
        excludedDoor = excludedDoor(randi([1,2])); % pick one of the empty doors randomly
    end
    
    % find the chance of winning by no door switching
    if myChoice==doorWithCar
        TotalWinsWithoutSwitch = TotalWinsWithoutSwitch + 1;
    end
    
    % find the chance of winning by door switching
    myChoice = AllDoors(AllDoors~=myChoice & AllDoors~=excludedDoor); % update my choice to the alternative door after host's exclusion of one of the doors.
    if myChoice==doorWithCar
        TotalWinsWithSwitch = TotalWinsWithSwitch + 1;
    end
    
    averageNumberOfWinsWithoutSwitch(iExperiment) = TotalWinsWithoutSwitch / iExperiment;
    averageNumberOfWinsWithSwitch(iExperiment) = TotalWinsWithSwitch / iExperiment;

end

% Now draw and export figures

figure('visible','on'); hold on; box on;
plot( averageNumberOfWinsWithoutSwitch ...
    , 'color', 'blue' ...
    , 'linewidth', 1.5 ...
    );
plot( averageNumberOfWinsWithSwitch ...
    , 'color', 'red' ...
    , 'linewidth', 1.5 ...
    );
axis([0,nExperiments,0,1]);
line1 = line( [1 nExperiments], [0.66666666 0.66666666] ...
            , 'LineStyle', '--' ...
            , 'color','green' ...
            , 'LineWidth',1.5 ...
            );
line2 = line( [1 nExperiments], [0.33333333 0.33333333] ...
            , 'LineStyle','--' ...
            , 'color','green' ...
            , 'LineWidth',1.5 ...
            );
xlabel('Experiment Repeat Number');
ylabel('Average Probability of Winning');

% add legend to figure
legend( {'No Door Switch','With Door Switch','theoretical prediction','theoretical prediction'} ...
      , 'Location', 'southeast' ...
      );
uistack(line1,'bottom');
uistack(line2,'bottom');

set(gca,'xscale','log');

set( gca ...
   , 'XMinorTick','on','YMinorTick','on' ...
   , 'FontSize', 13 ...
   , 'LineWidth', 1.25 ...
   );

saveas(gcf,'MontyGameResult.png');

%close(gcf);
Python

Here is an example Python implementation, monteGame.py,

#!/usr/bin/env python

import numpy as np
import matplotlib.pyplot as plt

ngames = 100000

chance_of_success_by_switching = np.zeros(ngames,dtype=np.double)
chance_of_success_by_not_switching = np.zeros(ngames,dtype=np.double)

options = [1,2,3]

for i in range(ngames):
    door_host_will_open = options.copy()
    door_with_car = np.random.randint(1,4)
    door_host_will_open.remove(door_with_car)
    door_i_choose = np.random.randint(1,4)
    if (door_i_choose!=door_with_car): door_host_will_open.remove(door_i_choose)
    
    door_to_switch = options.copy()
    door_to_switch.remove(door_i_choose)
    door_to_switch.remove(door_host_will_open[0])
    
    if i==0:
        if (door_to_switch[0]==door_with_car):
            chance_of_success_by_switching[i] = 1.0
        else:
            chance_of_success_by_not_switching[i] = 1.0
    else:
        if (door_to_switch[0]==door_with_car):
            chance_of_success_by_switching[i] = (chance_of_success_by_switching[i-1]*np.double(i) + 1.0)/np.double(i+1)
            chance_of_success_by_not_switching[i] = chance_of_success_by_not_switching[i-1]*np.double(i)/np.double(i+1)
        else:
            chance_of_success_by_not_switching[i] = (chance_of_success_by_not_switching[i-1]*np.double(i) + 1.0)/np.double(i+1)
            chance_of_success_by_switching[i] = chance_of_success_by_switching[i-1]*np.double(i)/np.double(i+1)

trials = np.linspace(1,ngames,ngames)
plt.semilogx(trials, chance_of_success_by_switching, 'r-')
plt.hold('on')
plt.semilogx(trials, chance_of_success_by_not_switching, 'b-')
plt.legend(['switch choice' , 'not switch choice'])
plt.xlabel('Game Number')
plt.ylabel('Average probability of winnig')
plt.title('{} Monte Hall Games: On average, you win if your switch!'.format(ngames))
plt.savefig('MontyGameResult.png')
plt.show()

Comments