/****************************************************************
* Chapter 4 Example 4.12, 4.13 in Johnson and Wichern 5th ed    *
* Assessing Normality assumption: Univariate and Bivariate      * 
*****************************************************************/
DM "output;clear;log;clear";
Options pagesize=45 linesize=80 PageNo=1 NoDate;
OPTIONS FORMCHAR="|----|+|---+=|-/\<>*"; 
	*To make your document display correctly without SAS Monospace font;
Title1 "Example 4.12 & 4.13";
ODS HTML body=    "norm-body.html" 
         contents="norm-contents.html" 
         frame=   "norm-frame.html" 
         page=    "norm-page.html"
         headtext="<title>Assessing Normality Assumption</title>"
         anchor="norm";
	* HTML files are to be saved in, for example, C:\Program Files\SAS\SAS 9.1;
title "Data in Prob 1.4";
Data p1_4;
 obs+1;
 Input X1-X3;
 Label 	X1="Sales"
		X2="Profits"
  		X3="Assets";
DataLines;
126974  4224  173297
96933  3835  160893
86656  3510  83219
63438  3758  77734
55264  3939  128344
50976  1809  39080
39069  2946  38528
36156  359  51038
35209  2480  34715
32416  2413  25636
;

Proc Print DATA=p1_4; Id Obs; run;
title "Array of Descriptive Statistics for X1 & X2";
Proc IML;			
	Use p1_4;
	Read ALL VAR{ X1 X2 } into X; N=NROW(X); 
	One=SHAPE(1,N,1);
	MeanVec=(One`*X)/N;			*Calculate Mean Vector;
	M=REPEAT(MeanVec,N,1); 		
	Sigma=(X-M)`*(X-M)/(N-1);	*Calculate Cov Matrix;
	Dis=VECDIAG((X-M)*INV(Sigma)*(X-M)`); *Calculate stat distance;
	print MeanVec; 
	print Sigma; 
	print Dis;

	CREATE D2 FROM Dis;
	APPEND FROM Dis;
run;
Title "X1 and X2 with distance";
Data pp1_4; merge p1_4 D2(RENAME=(COL1=D2)); 
	 Drop X3;
Proc Print data=pp1_4 NOOBS; 
run;
Proc Sort data=pp1_4 out=sort_d2; 
	 by D2;
run;
Title "Distances, normal quantiles, and chi-sq quantiles";
Data norm_p; set sort_d2; 
	NN=(_N_-.5)/10;
	Z=PROBIT(NN);		*Calculate Normal quantiles;
	Xsq=CINV(NN,2);		*Calculate Chi-sq quantiles;
	drop obs NN X1 X2;
Proc print data=norm_p;run;
Proc GPLOT data=norm_p;
	title 'Normal Probability Plot';
	Plot D2*Z / HAXIS=AXIS1 VAXIS=AXIS2 VREF=0;
		AXIS1 Label=('Standard Normal') ORDER=-2 to 2 by 1;
		AXIS2 Label=('Squared Distance');
		Footnote J=L 'Data Prob. 1.4'; run;
Proc GPLOT data=norm_p;
	title 'Chi-sq Probability Plot';
	Plot D2*Xsq / HAXIS=AXIS1 VAXIS=AXIS2 VREF=0;
		AXIS1 Label=('Chi-square') ORDER=0 to 7 by 1;
		AXIS2 Label=('Squared Distance');
		Footnote J=L 'Data Prob. 1.4'; 
run;

ODS HTML Close;