% LMcs130302_4th.m
% Larsen & Marx (2006) Introduction to Mathematical Statistics, 4th edition
% page 791. A case study solved by the paired t test
% Written by Eugene.Gallagher@umb.edu 11/14/10 Revised 11/16/10, 3/1/11
% X(:,1) Alamo rent-a-car, X(:,2)Avis rent-a-car
X=[48.99 51.99
    49.99 55.99
    42.99 47
    34.99 42.99
    42.99 44.95
    33.99 38.99
    59 69
    42.89 50.99
    47.99 49.99
    47.99 53.99
    35.99 42.99
    44.99 44.99]
D=(X(:,1)-X(:,2));    
[H,P,CI,STATS] = TTEST(X(:,1),X(:,2),0.05,'both');
fprintf(...
 'The mean ''Avis-Alamo'' difference in price was $%5.2f.\n',...
 mean(D));
fprintf(...
 'The sd for the difference in price was $%5.2f.\n',...
 STATS.sd);
fprintf(...
   'The CI for the difference in price was [$%5.2f $%5.2f]\n',...
 CI);
fprintf('The paired t statistic was %6.4f with %2.0f df\n',...
    STATS.tstat,STATS.df)
fprintf('The paired t test 2-tailed p=%6.4f\n',P);
[p,h,stats] = signtest(D,0,0.05,'method','exact');
fprintf('The sign test exact p=%6.4f\n',p);
[p,h,stats] = signtest(D,0,'method','approximate');
fprintf('The sign test approximate p=%6.4f\n',p);
[P,H,STATS] = signrank(X(:,1),X(:,2),'alpha',0.05,'method','exact');
fprintf('The signed rank statistic is %4.1f\n',STATS.signedrank);
fprintf('The sign rank test exact p=%6.4f\n',P);
[P,H,STATS] = signrank(X(:,1),X(:,2),'alpha',0.05,'method','approximate');
fprintf('The signed rank z statistic is %6.4f\n',STATS.zval);
fprintf('The sign rank test approximate p=%6.4f\n',P);
% Plot histogram and check for symmetry
binsize=2;  % Needed in order to properly scale the normal pdf
edges=-12:binsize:4;
% hist(D);
[N,BIN] = histc(D,edges);
bar(edges,N,'histc')
axis([-12.1 4.1 0 4.1])
set(get(gca,'Children'),'FaceColor',[.8 .8 1]);
xlabel('''Avis-Alamo'' Difference in Price','FontSize',20);
ylabel('Number of Cases','FontSize',20);
ax1=gca;
set(ax1,'Ytick',[0:4],'Xtick',[-12:2:4],'XTickLabel',...
    {'$12','$10','$8','$6','$4','$2','$0','$2','$4'})
figure(gcf);pause
% Superimpose the normal probability pdf on a histogram of differences.
% The normal probability equation is provided on p. 293
% This is for mean 0, and unit standard
% deviation. The more general equation (Legendre & Legendre, 1998 p. 147) is:
% f(y_j)=1/(sqrt(2*pi)*sigma_j)*exp(-1/2*((y_j-mu_j)/sigma_j)^2);
n=length(D);
mu_j=mean(D);
sigma_j=std(D); % sigmaj is the standard deviation; = 1 after Z transform
y_j=-12.1:0.1:4.1;
fy_j=1/(sqrt(2*pi)*sigma_j)*exp(-1/2*((y_j-mu_j)./sigma_j).^2);
fy_j=n*binsize*fy_j;  % will properly scale the height of the pdf
% 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',2);
set(h1,'linestyle','--','color','r','linewidth',2)
s=sprintf('Case Study 13.3.2, %1.0f samples untransformed',n);
title(s,'FontSize',22)
figure(gcf);pause
