%
% Written by:
% -- 
% John L. Weatherwax                2008-02-20
% 
% email: wax@alum.mit.edu
% 
% Please send comments and especially bug reports to the
% above email address.
% 
%-----

clear all; 
close all; 
clc; 

addpath('../../Code/CSTool');
load whisky;      % <- the total population 
fprintf('the population median for for the (PRIVATE,STATE) stores is given by=(%10.5f,%10.5f)\n', median(private), median(state)); 
nPop  = size(private,1)+size(state,1);
nPriv = size(private,1); indsPriv = 1:nPriv;
nSta  = size(state,1);   indsSta  = (nPriv+1):(nPriv+nSta); 

% the number of bootstrapped samples: 
B = 1000; 

alpha = 0.1; 

ssArray = [ 5,10,15,20,30,40 ]; % <- the size of the sample to consider ... 
for ssi=1:length(ssArray),
  ss   = ssArray(ssi);
  fprintf('with a sample size of=%10d the alpha=%10.5f\n',ss, alpha); 
  inds = unidrnd(nPop,1,ss); % indices from BOTH private and state populations ... 
  
  % get a samples from the private stores: 
  s_private = private( intersect(inds,indsPriv) ); 
  m_private = median( s_private );                 % <- the sample median for the private stores
  bootstat  = bootstrp( B, @median, s_private ); 
  p_ci      = quantile( bootstat, [ 0.5*alpha, 1-0.5*alpha ] ); 
  fprintf( 'The CI for the PRIVATE stores is given by (%10.5f,%10.5f)\n', p_ci(1), p_ci(2) );
  
  s_state  = state( intersect(inds,indsSta)-nPriv ); 
  m_state  = median( s_state );                     % <- the sample median for the private stores
  bootstat = bootstrp( B, @median, s_state ); 
  p_ci = quantile( bootstat, [ 0.5*alpha, 1-0.5*alpha ] ); 
  fprintf( 'The CI for the STATE stores is given by   (%10.5f,%10.5f)\n', p_ci(1), p_ci(2) );
  
end