% Solution to problem 1.11.2 in Yates and Goodman

% Number of simulations.  Here I am going to do 10000.  See below for why.
N = 10000;

% VERSION 1
% =========

% Allocate space for C an H
% Note that the vectors are initialized to zero so for H we need only
% record when a heads turns up
C1 = zeros(N,1);
H1 = zeros(N,1);
for n = 1:N
  
  % First randomly determine which coin is picked
  r = rand(1);
  if (r < 0.5)
    C1(n) = 1;
  else
    C1(n) = 2;
  end
  
  % Now determine whether we get a heads or tails from the chosen coin
  r = rand(1);
  if (C1(n) == 1)
    if (r < 0.75)
      H1(n) = 1;
    end
  else
    if (r < 0.5)
      H1(n) = 1;
    end
  end
  
end

% VERSION 2
% =========
% Same as version 1, but the logic for the head decision is a bit more
% terse

% Allocate space for C an H
% Note that the vectors are initialized to zero so for H we need only
% record when a heads turns up
C2 = zeros(N,1);
H2 = zeros(N,1);
for n = 1:N
  
  % First randomly determine which coin is picked
  r = rand(1);
  if (r < 0.5)
    C2(n) = 1;
  else
    C2(n) = 2;
  end
  
  % Now determine whether we get a heads or tails from the chosen coin
  r = rand(1);
  if ( ((C2(n) == 1) & (r<0.75)) | ((C2(n) == 2) & (r<0.5)) )
    H2(n) = 1;
  end
  
end

% VERSION 3
% =========
% Here we use the sample space at the end of the tree in the exmaple.  For
% each trial we generate a single random number between 0 and 1.  If that
% number is between 0 and 3/8, the result is C1 and Heads.  I the number is
% between 3/8 and 4/8, the outcome we record is C1 and Tails.  If the
% number is between 4/8 and 6/8, the outime is C2 and Heads.  Otherwise the
% result is C2 and Tails.  So each possible element of the sample spac is
% allocated a range of valies for the random number equal to the
% probability computed at the end of the tree.

C3 = zeros(N,1);
H3 = zeros(N,1);
for n = 1:N
  r = rand(1);
  if ( r < 3/8)
    C3(n) = 1;
    H3(n) = 1;
  elseif (r < 4/8)
    C3(n) = 1;
  elseif (r < 6/8)
    C3(n) = 2;
    H3(n) = 1;
  else
    C3(n) = 2;
  end
end


% To check each appproach, let's simulate not 50 but 10000 runs and see if
% we get the right fraction of heads (5/8=0.625 of the time) and percent
% time we choose C1 (half the time)

% Fractions of heads
heads = sum([H1 H2 H3])/N
% Fractions of C1
coin1= sum([C1==1 C2==1 C3==1])/N

% The results I get are
% heads =
%     0.6258    0.6231    0.6224
% coin1 =
%     0.4977    0.4921    0.4935

% So, it seems as though we are pretty close to 0.625 and 0.5.  Each
% approach seems to be working.  Later in the term we will determine
% whether the diffferences in the above results from what we expect are in
% a precise statistical sense "significant."