# Robust Wald test for Poisson with Stata

I have a Poisson model which I use on a medical data set with 329 observations. In the regression I am particularly interested in the simple Wald statistic for a single coefficient: $$t = \frac{\widehat{\beta} - \beta}{s.e.(\widehat{\beta})}$$ I was told that simply using the robust option will not help and it was suggested that I use the percentile t-method for conducting the Wald test to refine the t-statistic. I know that I have to use the bootstrap command in Stata but I have not figured it out completely. The theoretical points about this type of refinement are clear but the problem is to implement it and then to know what the resulting $t$ tells me compared to the un-refined one. Any guidance on this would be much appreciated.

As far as I know there is no ready made command for your purpose in Stata but it does not seem to be necessary because it can be easily implemented by hand. If you run your regression

poisson y x, vce(robust)


Create a local which holds the coefficient of your variable of interest (call it "bx", for instance) and then use the bootstrap command on the test statistic you posted.

bootstrap twald = ((_b[x] - bx')/_se[x]), reps(800) nodots: poisson y x, vce(robust)


If you need your results to be replicable, set a seed first. Once you have the result you can compare the bootstrapped standard error of "twald" with the standard deviation of the standard normal which is 1 (if this is what you meant by comparing the refined t with the un-refined t).

The Robust Wald test in Stata is achieved by using the option robust after the glm or poisson command in Stata. This is not the same as bootstrapped standard error estimates, which are another type of "robust" estimate for standard errors. The nominal "robust standard errors" estimates you receive by using the robust` command are the Huber White standard error estimates.

Huber White estimators are not resampling based, but use the empirical information to estimate $\widehat{\mbox{se}(\hat{\beta})}$. They account for unspecified sources of correlation in the data, heteroscedasticity, and model mispecification. The parameter estimates from such model coefficients are interpreted as "population averaged" which makes them useful for non-causal epidemiologic studies of association.