options linesize=76 pagesize=76;

data use;
infile cards;
input ask   sell time lot rms bdrms rec;

/*
ask: asking price; sell:selling price; time:days on sale;
lot:sq.ft. of lot; rms: nb of rms; bdms: nb of bedrooms; rec: has recroom?
*/
cards;
44900 42000 106 5850 6 3 0
41000 38500 39 4000 4 2 0
52900 49500 36 3060 6 3 0
62500 60500 83 6650 6 3 1
65000 61000 69 6360 5 2 0
69900 66000 44 4160 7 3 1
72900 66000 39 3880 7 3 0
72900 69000 12 4160 5 3 0
85900 83800 14 4800 5 3 1
93500 88500 80 5500 7 3 1
99900 90000 91 7200 8 3 0
31900 30500 74 3000 5 2 0
29900 27000 25 1700 6 3 0
37900 36000 7 2880 6 3 0
39900 37000 68 3600 4 2 0
39900 37900 11 3185 6 2 0
42900 40500 65 3300 6 3 0
54900 40750 78 5200 8 4 0
49900 45000 64 3450 5 1 0
49900 45000 49 3986 6 2 1
46900 48500 37 4785 6 3 1
71900 65900 24 4510 8 4 0
39900 37900 44 4000 6 3 0
41900 38000 42 3934 4 2 0
42900 42000 111 4960 4 2 0
43900 42300 68 3000 5 2 0
45900 43500 30 3800 5 2 0
46000 44000 58 4960 4 2 0
49900 44500 16 3000 5 3 0
49000 44900 33 4500 6 3 0
52900 45000 37 3500 6 2 0
57500 48000 119 3500 7 4 0
119900 113750 11 6000 9 3 1
124500 120000 79 7000 8 3 0
74500 70000 79 12900 6 3 0
75900 71000 77 7686 6 3 1
89500 82000 50 5000 7 3 0
85900 82000 63 5800 8 3 0
85900 82500 68 6000 7 3 0
87500 83000 49 4800 7 3 0
86900 84000 59 6500 7 3 0
88200 85000 22 7320 7 4 0
87900 85000 91 6525 7 3 0
94900 91500 37 4800 8 3 1
96900 94000 91 6000 7 3 0
109900 103000 174 6000 8 3 1
109900 105000 39 6000 8 3 1
109900 105000 21 6000 9 3 0
;
data final; set use; 
obs=_n_;
variance=1/ask;

title1 "A simple procedure to correct for hetero. that uses only the asking ";
title2 "price to capture hetero. effects"; 

proc reg data=final;
model sell= ask time lot rms bdrms rec/ acov dw spec;
weight variance;

data r;
title1 "The Feasible Aitken Procedure"; 
title2 "First collect the residuals in the output data set err";
title3 "and call the residuals resid."; 

proc reg data=final;
model sell= ask time lot rms bdrms rec/ acov dw spec;
output out=err r=resid;

data e2; set err;
e2=resid**2;

title2 "Now use the squared residuals e2 in a second regression to allow";
title3 "several variables at once to capture hetero. effects."; 
title4 "Notice that this regr. can be used to test for hetero. effects."; 

proc reg data=e2;
model e2 = ask time rms bdrms ;
output out=var p=pred_var;

data feasibl; merge final var;
variance=1/pred_var;

title2 "Now use the predicted squared residuals pred_var to transform the";
title3 "original regression to get rid of hetero.";
title4 " "; 

proc reg data=feasibl;
model sell= ask time lot rms bdrms rec  / acov dw spec;
weight variance;

/* The proc means below collects the average of the squared errors which will 
be used in the LM test */

proc means data=err noprint; var resid; output out=sigma var=s2;

data lm; set e2; if _n_=1 then set sigma;
lm=e2/s2;
ask2=ask**2;    at=ask*time; al=ask*lot;  ar=ask*rms;   ab=ask*bdrms;ac=ask*rec;
time2=time**2;  tl=time*lot; tr=time*rms; tb=time*bdrms;tc=time*rec;
lot2=lot**2;    lr=lot*rms;  lb=lot*bdrms;lc=lot*rec;
rms2=rms**2;    rb=rms*bdrms;rc=rms*rec;
bdrms2=bdrms**2;bc=bdrms*rec;

title1 "Some LM tests for Heterosckedasticity";
title2 "The Breusch-Pagan Test";
title3 " "; 
proc reg data=lm;
model lm = ask time lot rms bdrms rec
           ask2 time2 lot2 rms2 bdrms2
           at al ar ab ac 
              tl tr tb tc 
                 lr lb lc 
                    rb rc 
                       bc; 

run;
quit;