function [bootStat,repStat]=bootExample();

%Have STAT 5114 code this up themselves (second week of class)


N=50;  %data samples
R=1000; %replicates from population
%Note: Stat of interest is hard coded in two places.

B=500; %bootstrap samples 
mixtureWeight1=1/2;
mixtureWeight2=1/3;
mixtureWeight3=1/6;

%Sample from population
for(r=1:R)
    for(i=1:N)    
        u=rand;
        if(u<mixtureWeight1)
            X(i)=randn+10; 
        elseif(u<mixtureWeight1+mixtureWeight2)
            X(i)=randn-10;
        else
            X(i)=randn+3;
        end
    end
    figure(1)
    hist(X,100)
    title('Hist from population')
    pause(.001)
    %repStat(r) = mean(X);  %true iid replicated
    repStat(r) = var(X);
end
figure(2)
subplot(1,2,1)
hist(repStat)
title('stat replicated from pop')

%Note: we are recycling only the final X-vector here! 1 replicate
%-------- bootstrap-----------
for(b=1:B)
    XBoot = randsample(X,N,true);  %samples N times w/replacement
    bootStat(b) = mean(XBoot); %statistic of interest, with bootstrap replication
    bootStat(b) = var(XBoot);
end
figure(2)
subplot(1,2,2)
hist(bootStat)
title('bootstrap distribution')

figure(3)
qqplot(bootStat,repStat)
xlabel('BootStat Quantiles')
ylabel('Replicate Quantiles')
title('QQ Plot')
