/*
data from article <<Changes in Fatal and nonfatal
                    crash rates on a toll highway>>

      by Doege T and Levy PS. Amer J Epid 103 236-241, 1976.
*/

DATA sasuser.tolldata;

INPUT

year           /* 1968 ... 1974                             */

d_100mvm       /* d[enominator]...100 million vehicle miles */

c_no_inj       /* no. of crashes without injury             */
c_nfatal       /* no. of crashes with nonfatal injury       */
c_fatal        /* no. of crashes with fatality              */
c_total        /* no. of crashes ... total                  */

;

d_vm000 = d_100mvm * 100000; /* d[enominator]...1000 vehicle miles   */
                             /* need as integer if use as Binomial   */
                             /* e.g.  8.90937 x 10^8 = 890937 x 10^3 */

ln_d000 = log(d_vm000);      /* for use as offset in log_rate model  */


year_68 = year - 1968;  /* year zero = 1968            */
year_71 = year - 1971;  /* centered over middle year   */
year_74 = year - 1974;  /* year zero = 1974            */

I_change = 0; if year = 1974 then I_change = 1;

                       /* indicator of change in speed limit */

d_year68 = d_100mvm*year_68;      /* for use in additive-rate model  */

r_no_inj  = c_no_inj / d_100mvm;
r_nfatal  = c_nfatal / d_100mvm;
r_fatal   = c_fatal  / d_100mvm;
r_total   = c_total  / d_100mvm;

LINES;
1968   8.90937  627  290  17  934
1969   9.70208  878  353  19 1250
1970  10.17554  840  334  16 1190
1971  11.12855  951  361  16 1328
1972  11.09100 1487  472  19 2036
1973  12.40563 1252  394  12 1658
1974  12.69948  931  318   5 1254
;

title 'all crashes .. null model via proc logistic ';
run;

proc logistic;
  model  c_total / d_vm000 = ;
run;

title 'all crashes.. null logistic model via proc genmod ';
run;

proc genmod;
  model  c_total / d_vm000 = /
                link = logit dist = binomial obstats;
run;

title 'all crashes.. null logistic proc genmod overdispersion';
run;

proc genmod;
  model  c_total / d_vm000 = /
            link = logit dist = binomial dscale;
run;

title 'all crashes.. null logistic proc genmod  sandwich se';
run;

proc genmod;
  class year;
  model  c_total / d_vm000 = / link = logit dist = binomial;
  repeated subject = year / modelse;
run;



title 'all crashes.. null Poisson log-rate model via proc genmod ';
run;

proc genmod;
  model  c_total = / link = log  dist = poisson offset = ln_d000;
run;

title 'all crashes BEFORE THE CHANGE .. null Poisson log-rate model ';
run;

proc genmod; where (year < 1974);
  model  c_total = / link = log  dist = poisson offset = ln_d000;
run;

title 'all crashes BEFORE THE CHANGE .. null Logit model ';
run;

proc genmod; where (year < 1974);
  model  c_total / d_vm000= / link = logit  dist = binomial ;
run;

title 'all crashes BEFORE THE CHANGE .. non-null Poisson log-rate model ';
run;

proc genmod; where (year < 1974);
  model  c_total = year_68 / link = log  dist = poisson offset = ln_d000;
run;

title 'all crashes BEFORE THE CHANGE .. non-null Logit model ';
run;
proc genmod; where (year < 1974);
  model  c_total / d_vm000= year_68 / link = logit  dist = binomial ;
run;

title 'all crashes BEFORE THE CHANGE .. additive rates Poisson model ';
run;
proc genmod; where (year < 1974);
  model  c_total  = d_100mvm / link = identity  dist = poisson noint ;
run;

title 'all crashes BEFORE THE CHANGE .. additive rates Poisson model ';
run;
proc genmod; where (year < 1974);
  model  c_total  = d_100mvm d_year68 / link = identity  dist = poisson noint ;
run;