a=100;  %interval length
N=2000; %# of samples 
replicates=1000;   %# of times to repeat the experiment

for(i=1:replicates)
    x=rand(1,N).*a;          %Generate the data
    
    r(i)=max(x)-min(x);      %Range
    v(i)=(min(x)+max(x))/2;  %Midrange
    xBar(i)=mean(x);         %Mean

    weights=linspace(1/N,1,N);
    quant_ls(i)=weights*sort(x)' / (weights*weights');           %Least squares on the quantiles
    
    bayes_estimator(i)= ((N+1)/N)*max(x);%((N+2)/(N+1))*max(x);  %various Bayes Estimators
end
figure(1)
clf
subplot(1,5,1); hist(r);
title('$Range$','interpreter','latex')
subplot(1,5,2); hist(2.*v);
title('$2\times Midrange$','interpreter','latex')
subplot(1,5,3); hist(2.*xBar);
title('$2\times \bar x$','interpreter','latex')
subplot(1,5,4); hist(quant_ls);
title('$Least Squares Quantile Based regression$','interpreter','latex')
subplot(1,5,5); hist(bayes_estimator);
title('$Bayes Estimator$','interpreter','latex')

figure(2)
plot(r,v,'o')
xlabel('Range','fontsize',18)
ylabel('Midrange','fontsize',18)



fprintf('MSE(range): %11.2f\n',(mean(r) - a)^2 + var(r));
fprintf('MSE(2*midrange): %6.2f\n',(mean(2*v) - a)^2 + var(2*v));
fprintf('MSE(2*mean): %6.2f\n',(mean(2*xBar) - a)^2 + var(2*xBar));
fprintf('MSE(Quantile Least Squares): %6.2f\n',(mean(quant_ls) - a)^2 + var(quant_ls));
fprintf('MSE(Bayes Estimator): %6.2f\n',(mean(bayes_estimator) - a)^2 + var(bayes_estimator));