% LMex060202_4th.m
% Larsen & Marx (2006, p. 437).
% Introduction to Mathematical Statistics, 4th edition
% Normal probability pdf
% Written 11/3/10 by Eugene.Gallagher@umb.edu
% written for EEOS601
Ybar=502;   % The sample average
n=86;       % the sample size
muo=494;    % This mean is assumed to be a known parameter
sigmao=124; % This standard deviation is assumed to be a known parameter

z=(Ybar-muo)/(sigmao/sqrt(n));
fprintf('The z statistic is %6.4f\n',z)
alpha=0.05; twotailed=1; % change two-tailed to 0 for a 1-tailed test
if twotailed==1
    if z >= norminv(1-alpha/2) | z<=norminv(alpha/2)
        fprintf('Reject H_o at alpha=%4.2f level of significance, 2-tailed.\n',alpha);
    else
        fprintf(...
            'Fail to reject H_o at alpha=%4.2f level of significance, 2-tailed.\n',alpha)
    end
    if z<0
        pvalue=normcdf(z);
    else
        pvalue=1-normcdf(z);
    end
    fprintf('The two-sided p-value is %6.4f\n',pvalue*2)
    % This is also available as a statistics toolbox m.file, need to use
    % repmat to get n observations
    [H,P,CI,ZVAL] = ztest(repmat(Ybar,n,1),muo,sigmao,alpha,'both');
    fprintf('The two-tailed p-value is %6.4f\n',P);
else
    % this is for a right-tailed test, rewrite for a left tailed test or
    % use Matlab function with 'left' instead of 'right'
    if z > norminv(1-alpha)
        fprintf('Reject H_o at alpha=%4.2f level of significance, 1-tailed.\n',alpha);
    else
        fprintf(...
            'Fail to reject H_o at alpha=%4.2f level of significance, 1-tailed\n',alpha)
    end
    if z<0
        pvalue=normcdf(z);
    else
        pvalue=1-normcdf(z);
    end
    fprintf('The one-sided p-value is %6.4f\n',pvalue)
    % This is also available as a statistics toolbox m.file, need to use
    % repmat to get n observations
    [H,P,CI,ZVAL] = ztest(repmat(Ybar,n,1),muo,sigmao,alpha,'right');
    fprintf('The one-tailed p-value is %6.4f\n',P);
end

% plot Figure 6.2.5 % Delete this section if you are using this program to
% solve numeric homework problems
zlow=-abs(ZVAL);  % Make this statement general.
zup=abs(ZVAL);
z=-3.05:.001:3.05;
fzz=normpdf(z);
plot(z,fzz,'linestyle','--','color','r','linewidth',3)
ylabel('f_z(z)','FontSize',20)
xlabel('z','FontSize',22)
axis([-3.05 3.05 0 0.405])
set(gca,'Ytick',[0 0.4],'FontSize',18)
set(gca,'Xtick',[-abs(ZVAL) 0 abs(ZVAL)],...
    'XtickLabel',{'-0.60', '0','0.60'},'FontSize',18)
ax=axis;
ax1=gca; % save the handle of the graph
title('Figure 6.2.6','FontSize',22)
hold on
fz=normpdf([zlow zup]);
plot([zlow zlow;zup zup]',[0 fz(1);0 fz(2)]','-k','linewidth',1) 
% Fill in the upper tail with fill
y2=zup:.001:ax(2);
fy2=normpdf(y2);
fymax=normpdf(ax(2));
fill([zup y2 ax(2) ax(2)],[0 fy2 fymax 0],[.8 .8 1])
% Fill in the lower tail with fill
y3=ax(1):.001:zlow;
fy3=normpdf(y3);
fymin=normpdf(ax(1));
fill([ax(1) ax(1) y3 zlow],[0 fymin fy3 0],[.8 .8 1])
t=sprintf('Area = 0.2748');
text(ax(1)+.3,.3,t,'Color','b','FontSize',20);
text(1.3,.3,t,'Color','b','FontSize',20);
t1=sprintf('P-value = 0.2748+0.2748\n            = 0.5486');
text(0.6,0.38,t1,'Color','b','FontSize',20)
figure(gcf)
hold off