/* This program is intended to illustrate the modeling of intervention
   data ala Box and Tiao (1975)  */

/*  This program is taken from the textbook "Forecasting Examples" (SAS Institute, 1996)
    pp. 11-12 and Example 6 pp. 101 - 110.  */

data drivers;
   format date monyy.;
   input date:monyy5. injuries @@;
   beltlaw = (date ge '01jan83'd);
   title 'Driver Casualties Data';
   title2 'monthly totals';
   datalines;
jan80  1665  feb80  1361  mar80  1506  apr80  1360
may80  1453  jun80  1522  jul80  1460  aug80  1552
sep80  1548  oct80  1827  nov80  1737  dec80  1941
jan81  1474  feb81  1458  mar81  1542  apr81  1404
may81  1522  jun81  1385  jul81  1641  aug81  1510
sep81  1681  oct81  1938  nov81  1868  dec81  1726
jan82  1456  feb82  1445  mar82  1456  apr82  1365
may82  1487  jun82  1558  jul82  1488  aug82  1684
sep82  1594  oct82  1850  nov82  1998  dec82  2079
jan83  1494  feb83  1057  mar83  1218  apr83  1168
may83  1236  jun83  1076  jul83  1174  aug83  1139
sep83  1427  oct83  1487  nov83  1483  dec83  1513
jan84  1357  feb84  1165  mar84  1282  apr84  1110
may84  1297  jun84  1185  jul84  1222  aug84  1284
sep84  1444  oct84  1575  nov84  1737  dec84  1763
jan85  .     feb85  .     mar85  .     apr85  .
may85  .     jun85  .     jul85  .     aug85  .
sep85  .     oct85  .     nov85  .     dec85  .
;

goptions cback=white colors=(black) border reset=(axis symbol);

axis1 offset=(1 cm) label=('Year')
      order=('01jan80'd to '01jan85'd by year);
axis2 label=(angle=90 'Casualties')
      order=(800 to 2200 by 200);

symbol1 i=join l=1 v=star;

proc gplot data=drivers;
   format date year4.;
   plot injuries*date=1 / haxis=axis1
                          vaxis=axis2
                          vminor=1
                          href='01jan83'd
                          lh=2;
run;


/* The following program is found on pp. 103-109 of   */
/* Forecasting Examples for Business and              */
/* Economics Using the SAS System                     */

/* pp. 103-104 */
%dftest(drivers,injuries,dlag=12)
%put &dfpvalue;

proc arima data=drivers;
   i var=injuries crosscorr=beltlaw noprint;
   e input=beltlaw plot;
run;

/* p. 106 */
   e input=beltlaw p=(1)(2 12);
run;

/* p. 108 */
   f lead=12 out=fore2
             id=date
             interval=month;
run;

/* p. 109 */
goptions cback=white colors=(black) border reset=(axis symbol);

axis1 offset=(1 cm) label=('Year')
      order=('01jan80'd to '01jan86'd by year);
axis2 label=(angle=90 'Casualties')
      order=(750 to 2250 by 500);

symbol1 i=join l=1;
symbol2 h=2 pct v=F;
symbol3 l=20 i=join;

proc gplot data=fore2;
   format date year4.;
   plot injuries*date=1
        forecast*date=2
        l95*date     =3
        u95*date     =3 / overlay
                          haxis=axis1
                          vaxis=axis2
                          vminor=4
                          href='01jan83'd
                          lh=2;
run;