/* SAS program to accompany the blog post
"Head-tail versus head-head: A counterintuitive property of coin tosses"
by Rick Wicklin, published 13APR2016
http://blogs.sas.com/content/iml/2016/04/09/counterintuitive-property-of-coin-tosses
*/
/*
Story from Wired magaxine:
http://www.wired.com/2016/03/mathematicians-discovered-prime-conspiracy/
"A lecture at Stanford by the mathematician Tadashi Tokieda, of the
University of Cambridge, in which he mentioned a counterintuitive
property of coin-tossing: If Alice tosses a coin until she sees a
head followed by a tail, and Bob tosses a coin until he sees two
heads in a row, then on average, Alice will require four tosses
while Bob will require six tosses (try this at home!), even though
head-tail and head-head have an equal chance of appearing after two
coin tosses."
*/
data CoinFlip(keep= player toss);
call streaminit(1234567);
label Player = "Stopping Condition";
do i = 1 to 100000;
do player = "HT", "HH"; /* first Alice tosses; then Bob */
first = rand("Bernoulli", 0.5); /* 0 is Tails; 1 is Heads */
toss = 1; done = 0;
do until (done);
second = rand("Bernoulli", 0.5);
toss = toss + 1; /* how many tosses until done? */
if player = "HT" then /* if Alice, ... */
done = (first=1 & second=0); /* stop when HT appears */
else /* if Bob, ... */
done = (first=1 & second=1); /* stop when HH appears */
first = second; /* remember for next toss */
end;
output;
end;
end;
run;
proc means data=CoinFlip MAXDEC=1 mean median Q3 P90 max;
class player;
var toss;
run;
/* comparative histograms */
ods graphics / height=500px width=600px;
title "Distribution of Coin Tosses until HT or HH Observed";
title2 "100,000 Games Simulated";
proc sgpanel data=CoinFlip;
panelby Player / rows=2 layout=rowlattice rowheaderpos=left;
histogram Toss / binwidth=1;
colaxis values=(2 to 30);
rowaxis grid;
run;
/* comparative ECDFs */
proc univariate data=CoinFlip;
label Player = "Stopping Condition";
class player;
cdfplot toss / overlay vref=50 75 90;
ods select cdfplot;
run;
/* simulate 100,000 games. Keep track of how many Alice and Bob win */
data Game(keep=winner);
call streaminit(1234567);
do i = 1 to 100000;
first = rand("Bernoulli", 0.5); /* 0 is Tails; 1 is Heads */
done = 0;
do until (done);
second = rand("Bernoulli", 0.5);
if (first=1 & second=0) then do;
done = 1;
winner = "HT"; /* stop when HT appears */
end;
if (first=1 & second=1) then do;
done = 1;
winner = "HH"; /* stop when HH appears */
end;
first = second; /* remember for next toss */
end;
output;
end;
run;
proc freq data=Game;
tables winner / nocum binomial(level='HT');
run;