-------------------------------------------------------
*****   Zero-Inflated Poisson Simulations         *****
-------------------------------------------------------

SETUP:  we simulated data sets of size N=150 observations 
        and fit poisson regression to obtain a estimate
        for the regression coefficient.

        To compute standard errors of the regression estimate
        we considered:
           1.  model-based assuming Poisson is correct
           2.  sandwich assuming var(y) = phi * mu
           3.  sandwich assuming var(y) = mu + alpha * mu^2
           4.  empirical sandwich variance estimator

        In addition, we computed the negative binomial MLE
        based on the NB-2 parameterization.  For this we
        considered standard errors based on:
           1.  model-based assuming NB-2
           2.  empirical sandwich variance estimator

RESULTS:

(a) For the Poisson regression estimator


Standard error estimates:  (mean of 1000 estimates)

       beta  mbeta   seT  seM1  seM2  seM3  seM4
Int  -0.511 -0.531 0.207 0.133 0.227 0.216 0.202
x1    2.288  2.293 0.496 0.253 0.429 0.563 0.490
x2    1.000  1.011 0.289 0.172 0.291 0.308 0.287

     mbeta = mean estimate of beta
     seT   = variance of beta-hat across 1000 simulations
     seM1  = mean s.e. estimate assuming Poisson
     seM2  = mean s.e. estimate assuming scaled variance (QL)
     seM3  = mean s.e. estimate assuming NB-2 variance using MoM for alpha
     seM4  = mean s.e. estimate using emprical variance

(b) For the Negative Binomial MLE (NB-2)

Standard error estimates:

> print( round( out, 3 ) )
       beta mbeta.nb seT.nb  seNB seM4.nb
Int  -0.511   -0.534  0.205 0.243   0.199
x1    2.288    2.301  0.492 0.663   0.491
x2    1.000    1.008  0.286 0.376   0.286


(c) Coverage properties for the POISSON model estimates given in (a):

Coverage:

       beta covM1 covM2 covM3 covM4
Int  -0.511 0.801 0.964 0.952 0.944
x1    2.288 0.688 0.914 0.967 0.947
x2    1.000 0.750 0.950 0.963 0.943