FOCUS AREAS

SAS/STAT Examples

Fitting Bayesian Zero-Inflated Poisson Regression Models with the MCMC Procedure


Contents | SAS Program | PDF

  
data catch;
   input gender $ age count @@;   
   if gender = 'F' then do;
      female = 1; male = 0;
      end;
   else do;
      female = 0; male = 1;
      end;
   obs = _N_;
datalines;
   F 54 18 M 37 0 F 48 12 M 27 0
   M 55 0 M 32 0 F 49 12 F 45 11
   M 39 0 F 34 1 F 50 0 M 52 4
   M 33 0 M 32 0 F 23 1 F 17 0
   F 44 5 M 44 0 F 26 0 F 30 0
   F 38 0 F 38 0 F 52 18 M 23 1
   F 23 0 M 32 0 F 33 3 M 26 0
   F 46 8 M 45 5 M 51 10 F 48 5
   F 31 2 F 25 1 M 22 0 M 41 0
   M 19 0 M 23 0 M 31 1 M 17 0
   F 21 0 F 44 7 M 28 0 M 47 3
   M 23 0 F 29 3 F 24 0 M 34 1
   F 19 0 F 35 2 M 39 0 M 43 6
;


ods graphics on;
proc mcmc data=catch seed=1181 nmc=100000 thin=10
   propcov=quanew monitor =(_parms_ Pearson);
   ods select Parameters PostSummaries PostIntervals tadpanel;
   parms beta0 0 beta1 0 beta2 0;
   prior beta: ~ normal(0,var=1000);
   mu = exp(beta0 + beta1*female*age + beta2*male*age);
   model count ~ poisson(mu);
   if obs = 1 then Pearson = 0;
   Pearson =  Pearson + ((count - mu)**2/mu);
run;     
ods graphics off;  

ods graphics on;
proc mcmc data=catch seed =1181 nmc=100000 thin=10 
   propcov=quanew monitor =(_parms_ Pearson);
   ods select Parameters PostSummaries PostIntervals tadpanel;
   parms beta0 0 beta1 0 beta2 0 eta .3;
   prior beta: ~ normal(0,var=1000);
   prior eta ~ uniform(0,1);
   mu=exp(beta0 + beta1*female*age + beta2*male*age);           
   llike=log(eta*(count eq 0) + (1-eta)*pdf("poisson",count,mu));
   model general(llike);
   if obs = 1 then Pearson = 0;
   mean = (1 - eta)*mu;
   sigma2 = (1 - eta)*mu*(1 + eta*mu);   
   Pearson =  Pearson + ((count - mean)**2/sigma2);
run;
ods graphics off;