The CNTSELECT Procedure
Example 9.2 Zero-Inflated Poisson Model with BAYES and PRIOR Statements
This example shows how to use the CNTSELECT procedure to estimate a zero-inflated Poisson model by using Bayesian methods, with user-specified priors. The following DATA step generates 1,000 replicates from the zero-inflated Poisson (ZIP) model. The model contains seven variables and three variables that correspond to the zero-inflated process.
data bayes_ex;
call streaminit(12345);
array vars x1-x7;
array zero_vars z1-z3;
array parms{7} (.3 .4 .2 .4 -.3 -.5 -.3);
array zero_parms{3} (-.6 .3 .2);
intercept=0.5;
group=1;
z_intercept=-1;
theta=0.5;
do i=1 to 1000;
sum_xb=0;
sum_gz=0;
if i>500 then do;
intercept=2;
group=2;
end;
do j=1 to 7;
vars[j]=rand('NORMAL',0,1);
sum_xb=sum_xb+parms[j]*vars[j];
end;
mu=exp(intercept+sum_xb);
y_p=rand('POISSON', mu);
do j=1 to 3;
zero_vars[j]=rand('NORMAL',0,1);
sum_gz = sum_gz+zero_parms[j]*zero_vars[j];
end;
z_gamma = z_intercept+sum_gz;
pzero = cdf('LOGISTIC',z_gamma);
cut=rand('UNIFORM');
if cut<pzero then y_p=0;
output;
end;
keep y_p group x1-x7 z1-z3;
run;
You can load the bayes_ex data set into your CAS session by naming your CAS engine libref in the first statement of the following DATA step:
data mycas.bayes_ex;
set bayes_ex;
run;
These statements assume that your CAS engine libref is named mycas, but you can substitute any appropriately defined libref.
The following statements estimate a zero-inflated Poisson model via Bayesian analysis. Note that the BAYES statement controls the settings for the MCMC sampler and determines which Bayesian output to produce. The PRIOR statements set the prior distributions for the parameters.
proc cntselect data=mycas.bayes_ex dist=zip;
class group;
model y_p=group x1-x7;
zeromodel y_p ~ z1-z3;
bayes seed = 81239 nmc = 10000 sampler = rwm(ntu = 50) priorsummary(shownames);
prior x1-x7 ~ normal(mean = 0, sd = 10);
prior group_1 ~ normal(mean = 0, sd = 10, upper = 0);
prior z1-z3 ~ normal(mean = 0, sd = 1);
prior Intercept ~ normal(mean = 0, sd = 10);
run;
Output 9.2.1 shows the results for the zero-inflated Poisson model. The "Prior Summary" table shows detailed information about prior distributions for the parameters in the model. The "Posterior Summaries" table contains two point estimates of the parameters from the posterior samples: the posterior mean and the posterior median (50th percentile). The "MCMC Diagnostic Summaries" table contains basic convergence diagnostics to check whether the Markov chain has converged and whether the sample size is sufficient.
Output 9.2.1: Bayesian Estimation of Zero-Inflated Poisson Model
| Prior Summary | ||||||||
|---|---|---|---|---|---|---|---|---|
| Parameter | Name | Prior | Bounds | Hyperparameters | ||||
| Lower | Upper | Name | Value | Name | Value | |||
| Intercept | Intercept | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| group 1 | group_1 | Truncated Normal | -Infty | 0 | Mean | 0 | Variance | 100 |
| x1 | x1 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| x2 | x2 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| x3 | x3 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| x4 | x4 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| x5 | x5 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| x6 | x6 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| x7 | x7 | Normal | -Infty | Infty | Mean | 0 | Variance | 100 |
| Inf_Intercept | Inf_Intercept | Normal | -Infty | Infty | Mean | 0 | Variance | 1000000 |
| Inf_z1 | Inf_z1 | Normal | -Infty | Infty | Mean | 0 | Variance | 1000000 |
| Inf_z2 | Inf_z2 | Normal | -Infty | Infty | Mean | 0 | Variance | 1000000 |
| Inf_z3 | Inf_z3 | Normal | -Infty | Infty | Mean | 0 | Variance | 1000000 |
| Posterior Summaries | ||||||||
|---|---|---|---|---|---|---|---|---|
| Parameter | N | Mean | Standard Deviation | Percentiles | ||||
| 2.5% | 25% | 50% | 75% | 97.5% | ||||
| Intercept | 10000 | 1.9939 | 0.0197 | 1.9554 | 1.9808 | 1.9939 | 2.0067 | 2.0350 |
| group 1 | 10000 | -1.4420 | 0.0367 | -1.5075 | -1.4696 | -1.4439 | -1.4159 | -1.3682 |
| x1 | 10000 | 0.2862 | 0.0155 | 0.2571 | 0.2756 | 0.2859 | 0.2967 | 0.3158 |
| x2 | 10000 | 0.4182 | 0.0146 | 0.3883 | 0.4082 | 0.4187 | 0.4285 | 0.4470 |
| x3 | 10000 | 0.2147 | 0.0162 | 0.1826 | 0.2040 | 0.2147 | 0.2265 | 0.2454 |
| x4 | 10000 | 0.4035 | 0.0141 | 0.3744 | 0.3946 | 0.4042 | 0.4131 | 0.4293 |
| x5 | 10000 | -0.2989 | 0.0185 | -0.3349 | -0.3123 | -0.2978 | -0.2853 | -0.2654 |
| x6 | 10000 | -0.5106 | 0.0137 | -0.5374 | -0.5210 | -0.5097 | -0.5009 | -0.4863 |
| x7 | 10000 | -0.3056 | 0.0147 | -0.3347 | -0.3157 | -0.3059 | -0.2954 | -0.2767 |
| Inf_Intercept | 10000 | -1.0511 | 0.0932 | -1.2210 | -1.1136 | -1.0517 | -0.9953 | -0.8465 |
| Inf_z1 | 10000 | -0.6225 | 0.0952 | -0.8018 | -0.6884 | -0.6240 | -0.5592 | -0.4266 |
| Inf_z2 | 10000 | 0.3959 | 0.0920 | 0.2228 | 0.3376 | 0.3900 | 0.4633 | 0.5777 |
| Inf_z3 | 10000 | 0.2824 | 0.0882 | 0.1176 | 0.2224 | 0.2783 | 0.3495 | 0.4428 |
| MCMC Diagnostic Summaries | ||||||
|---|---|---|---|---|---|---|
| Parameter | MCSE | MCSE/SD | ESS | Autocorrelation Time | ESS/N | |
| Intercept | * | 0.00285 | 0.1444 | 47.9410 | 208.6 | 0.00479 |
| group 1 | * | 0.00598 | 0.1629 | 37.7059 | 265.2 | 0.00377 |
| x1 | * | 0.00212 | 0.1374 | 52.9928 | 188.7 | 0.00530 |
| x2 | * | 0.00154 | 0.1053 | 90.1485 | 110.9 | 0.00901 |
| x3 | * | 0.00217 | 0.1334 | 56.1771 | 178.0 | 0.00562 |
| x4 | * | 0.00188 | 0.1334 | 56.1891 | 178.0 | 0.00562 |
| x5 | * | 0.00302 | 0.1628 | 37.7233 | 265.1 | 0.00377 |
| x6 | * | 0.00181 | 0.1319 | 57.4459 | 174.1 | 0.00574 |
| x7 | * | 0.00137 | 0.0936 | 114.1 | 87.6691 | 0.0114 |
| Inf_Intercept | 0.00940 | 0.1008 | 98.3253 | 101.7 | 0.00983 | |
| Inf_z1 | * | 0.0122 | 0.1286 | 60.4869 | 165.3 | 0.00605 |
| Inf_z2 | * | 0.00758 | 0.0824 | 147.3 | 67.8904 | 0.0147 |
| Inf_z3 | 0.00736 | 0.0835 | 143.3 | 69.7857 | 0.0143 | |
| *Autocorrelation Remains | ||||||