
close all force
clear all
 disp('The 1954 Polio Vaccine Trial')
  lw = 2; 
  set(0, 'DefaultAxesFontSize', 16);
  fs = 16;
  msize = 10;
  rand('seed',0)
%Example 1.4 (page 5): The 1954  Polio Vaccine Trial  

vaccine = [ones(1, 33), zeros(1, 200745-33)]';
placebo = [ones(1, 115) zeros(1,201229-115)]';
rrs=[];

B=5000 

for i = 1:B
    pvac  = sum(bootsample(vaccine))/200745;
    pplac = sum(bootsample(placebo))/201229;
    rr=pvac/pplac;
    rrs =[rrs rr];
end
rrstar025 = quantile(rrs,0.025); %
rrstar975 = quantile(rrs,0.975); % 
 [rrstar025  rrstar975]  %   0.1865    0.4134  (for B=50000)
 %
histn(rrs,0.1,0.01,0.55);
hold on
rrt = (sum(vaccine)/200745)/(sum(placebo)/201229);  % 0.2876
plot([rrt rrt],[0 10], 'r-','linewidth',3)
plot([rrstar025 rrstar025],[0 10], 'g-','linewidth',3)
plot([rrstar975 rrstar975],[0 10], 'g-','linewidth',3)
axis([0.1, 0.55,0, 9])
hold off
%print -depsc 'C:\STAT\Boot\Booteps\poliovaccine.eps'

 
