% LMFig070201_4th.m
% Use as a model LMFig040301_4th.m
% Page 471-472 in
% Larsen & Marx (2006) Introduction to Mathematical Statistics, 4th edition
% written by Eugene.Gallagher@umb.edu, 3/12/11 6:50, revised 3/12/11
% Generate 400 samples from a standard normal distribution, size 4
clf;cla;hold off
K=menu('Choose an Option','Figure 7.2.1a','Figure7.2.1b')
Binsize=.5;
n=4;N=500;
if K==1
    R1 = normrnd(0,1,n,N);
    Z1=(mean(R1)-0)./(1./sqrt(n));
    tl='Figure 7.2.1a';
else
    R2 = normrnd(0,1,n,N);
    Z1=(mean(R2)-0)./(std(R2)./sqrt(n));
    tl='Figure 7.2.1b';
end
edges=-4:Binsize:4;
[N1,BIN1] = histc(Z1,edges);
bar(edges,N1/N,'histc');
axis([-4.2 4.2 0 0.25]);
set(get(gca,'Children'),'FaceColor',[.8 .8 1])
ylabel('Density','FontSize',20);
ax1=gca;set(ax1,'xtick',-4:4,'FontSize',20);
set(ax1,'ytick',[0.1 0.2],'FontSize',20)
figure(gcf);pause
% Superimpose the normal probability pdf:
% The more general equation(Legendre & Legendre, 1998 p. 147, 
% or LM eq 7.1.1 is:
% f(y_j)=1/(sqrt(2*pi)*sigma_j)*exp(-1/2*((y_j-mu_j)/sigma_j)^2)
mu_j=0;
sigma_j=1; % sigmaj is the standard deviation
y_j=-4:0.01:4;
%fy_j=1/(sqrt(2*pi)*sigma_j)*exp(-1/2*((y_j-mu_j)./sigma_j).^2);
% or the identical function:
fy_j=normpdf(y_j,mu_j,sigma_j)*Binsize;
% fyj=1/(sqrt(2*pi)*sigmaj)*exp(-1/2*((y-muj)/sigmaj).^2); 
% Plot using ax1 handle, saved above,to save this graph
% on top of the previous graph.
h1=line(y_j,fy_j,'Color','r','Parent',ax1,'Linewidth',3);
set(h1,'linestyle','--','color','r','linewidth',3)
xlabel(''),ylabel('Density','FontSize',20)
title(tl,'FontSize',22)
figure(gcf);pause
hold off