options ls=72 ps=72;

/*------------------------------------------------------*/
/*                                                      */
/* Proc Genmod to fit population average (marginal)     */
/* model using GEE approach for the epileptic seizure   */
/* count data                                           */
/*                                                      */
/*------------------------------------------------------*/

data seizure; 
  infile "seize.dat";
  input id seize visit trt age;
  nobs=_n_;
  interval = 2;
  if visit=0 then interval=8;
  logtime = log(interval);
  assign = (visit>0);
run;


title "Model 1: overall effct of the treatment";
proc genmod data=seizure;
  class id;
  model seize = assign trt assign*trt 
     / dist=poisson link=log offset=logtime;
  repeated subject=id / type=exch;
run;


title "Model 2: take overdispersion into account";
proc genmod data=seizure;
  class id;
  model seize = assign trt assign*trt 
     / dist=poisson link=log offset=logtime scale=pearson aggregate=nobs;
  repeated subject=id / type=exch;
run;


title "Model 3: take randomization into account";
proc genmod data=seizure;
  class id;
  model seize = assign assign*trt 
     / dist=poisson link=log offset=logtime scale=pearson aggregate=nobs;
  repeated subject=id / type=exch;
run;


title "Model 4: adjusting for other covariates";
proc genmod data=seizure;
  class id;
  model seize = assign trt assign*trt age
     / dist=poisson link=log offset=logtime scale=pearson aggregate=nobs;
  repeated subject=id / type=exch;
run;