/* mvdeath.sas */
/* data for 49 states omitting DC and HI */
options ls=72 ps=60;
data mvdeath;
   infile 'mvdeath.dat';
   input  code   state $ drivers  rural    fuel  deaths;
   ldrivers = log(drivers);
   lrural = log(rural);
   lfuel = log(fuel);
   ldeaths = log(deaths);
proc plot data=mvdeath;
   plot deaths*drivers = '*';
   plot deaths*rural = '*';
   plot deaths*fuel = '*';
   plot drivers*rural = '*';
   plot drivers*fuel = '*';
   plot rural*fuel = '*';
title 'plots for untransformed variables';
proc plot data=mvdeath;
   plot ldeaths*ldrivers = '*';
   plot ldeaths*lrural = '*';
   plot ldeaths*lfuel = '*';
   plot ldrivers*lrural = '*';
   plot ldrivers*lfuel = '*';
   plot lrural*lfuel = '*';
title 'plots for log transformed variables';
/* fit models */
/* regression of deaths on drivers, rural, and fuel (full model) */
proc glm data=mvdeath;
   model deaths = drivers rural fuel  / ss1 ss2 clparm;
      output out=new1 predicted=yhat residual=resid
          lcl=clpl ucl=clpu lclm=clml uclm=clmu;
title 'regression of deaths on drivers, rural, and fuel (full model)';
proc print data=new1;
   var state deaths drivers rural fuel;
proc print data=new1;
   var state clml clmu clpl clpu yhat resid;
proc plot data=new1;
   plot resid*yhat / vref=0;
   plot resid*drivers / vref=0;
   plot resid*rural / vref=0;
   plot resid*fuel / vref=0;
title 'residual plots: full model untransformed';
/* regression of log(deaths) on log(drivers), log(rural), and log(fuel) (full model) */
proc glm data=mvdeath;
   model ldeaths = ldrivers lrural lfuel  / ss1 ss2 clparm;
      output out=new2 predicted=yhat residual=resid
          lcl=clpl ucl=clpu lclm=clml uclm=clmu;
title 'regression of ldeaths on ldrivers, lrural, and lfuel (full model)';
proc print data=new2;
   var state ldeaths ldrivers lrural lfuel;
proc print data=new2;
   var state clml clmu clpl clpu yhat resid;
proc plot data=new2;
   plot resid*yhat / vref=0;
   plot resid*ldrivers / vref=0;
   plot resid*lrural / vref=0;
   plot resid*lfuel / vref=0;
title 'residual plots: full model log transformed';
/* regression of log(deaths) on log(drivers) and log(rural) (reduced model) */
proc glm data=mvdeath;
   model ldeaths = ldrivers lrural  / ss1 ss2 clparm;
      output out=new3 predicted=yhat residual=resid
          lcl=clpl ucl=clpu lclm=clml uclm=clmu;
title 'regression of ldeaths on ldrivers and lrural (reduced model)';
proc print data=new3;
   var state ldeaths ldrivers lrural;
proc print data=new3;
   var state clml clmu clpl clpu yhat resid;
proc plot data=new3;
   plot resid*yhat / vref=0;
   plot resid*ldrivers / vref=0;
   plot resid*lrural / vref=0;
title 'residual plots: reduced model log transformed';
/* regression of log(deaths) on log(drivers) (one variable models) */
proc glm data=mvdeath;
   model ldeaths = ldrivers  / ss1 ss2 clparm;
      output out=new4 predicted=yhat residual=resid
          lcl=clpl ucl=clpu lclm=clml uclm=clmu;
title 'regression of ldeaths on ldrivers (one variable model)';
proc print data=new4;
   var state ldeaths ldrivers  clml clmu clpl clpu yhat resid;
proc plot data=new4;
   plot resid*ldrivers / vref=0;
   plot resid*lrural / vref=0;
title 'residual plots: one variable model log transformed';
proc glm data=mvdeath;
   model ldeaths = lrural  / ss1 ss2 clparm;
      output out=new5 predicted=yhat residual=resid
          lcl=clpl ucl=clpu lclm=clml uclm=clmu;
title 'regression of ldeaths on lrural (one variable model)';
proc print data=new5;
   var state ldeaths lrural  clml clmu clpl clpu yhat resid;
proc plot data=new5;
   plot resid*lrural / vref=0;
   plot resid*ldrivers / vref=0;
title 'residual plots: one variable model log transformed';
proc glm data=mvdeath;
   model ldeaths = lfuel  / ss1 ss2 clparm;
      output out=new6 predicted=yhat residual=resid
          lcl=clpl ucl=clpu lclm=clml uclm=clmu;
title 'regression of ldeaths on lfuel (one variable model)';
proc print data=new6;
   var state ldeaths lfuel  clml clmu clpl clpu yhat resid;
proc plot data=new6;
   plot resid*lfuel / vref=0;
   plot resid*lrural / vref=0;
   plot resid*ldrivers / vref=0;
title 'residual plots: one variable model log transformed';
/* residual analysis */
proc univariate data=new1 plot normal;
   var resid;
title 'residual analysis: model deaths = drivers rural fuel';
proc rank data=new1 normal=blom out=new1a;
   var resid;
   ranks nscore;
proc plot data=new1a;
   plot nscore*resid;
proc univariate data=new2 plot normal;
   var resid;
title 'residual analysis: model ldeaths = ldrivers lrural lfuel';
proc rank data=new2 normal=blom out=new2a;
   var resid;
   ranks nscore;
proc plot data=new2a;
   plot nscore*resid;
proc univariate data=new3 plot normal;
   var resid;
title 'residual analysis: model ldeaths = ldrivers lrural';
proc rank data=new3 normal=blom out=new3a;
   var resid;
   ranks nscore;
proc plot data=new3a;
   plot nscore*resid;
proc univariate data=new4 plot normal;
   var resid;
title 'residual analysis: model ldeaths = ldrivers';
proc rank data=new4 normal=blom out=new4a;
   var resid;
   ranks nscore;
proc plot data=new4a;
   plot nscore*resid;
proc univariate data=new5 plot normal;
   var resid;
title 'residual analysis: model ldeaths = lrural';
proc rank data=new5 normal=blom out=new5a;
   var resid;
   ranks nscore;
proc plot data=new5a;
   plot nscore*resid;
proc univariate data=new6 plot normal;
   var resid;
title 'residual analysis: model ldeaths = lfuel';
proc rank data=new6 normal=blom out=new6a;
   var resid;
   ranks nscore;
proc plot data=new6a;
   plot nscore*resid;
run;