Re: st: robust poisson regression vs. glm with log link

[Maarten Buis <maartenlbuis@gmail.com>]

On Tue, Aug 23, 2011 at 8:03 PM, Dimitriy V. Masterov wrote:
> William Gould has a new post on the NEC blog about how Poisson
> regression with vce(robust) is a better alternative to log-linear
> regression:
>
> http://blog.stata.com/2011/08/22/use-poisson-rather-than-regress-tell-a-friend/#disqus_thread
>
> I can't seem to comment there, so I hope I can post my questions here
> to get some insight. First, how does the robust Poisson model perform
> relative to glm with a log link function? Are there ways to get the
> elasticities from the poisson and/or glm coefficients?

You can get exactly the same estimates and standard errors by also
specifying the -family(poisson)- option. If you leave that option out
you are using a different variance function, in this case constant
variance. The logic comes from the quasi-likelihood tradition, which
has found that one can make a function that has very similar
properties to the likelihood function by just specifying how the mean
depends on the explanatory variables (the link-function) and how the
variance depends on that mean (the variance function). If all we want
to do is see how the mean changes when the explanatory variables
change, than there is an asymptotic argument that all that has to be
correctly specified is this mean function. In finite samples the
variance function may still matter, and that is why quite some
attention is typically spent on also checking the variance function
and why -glm- allows various variance functions (-family()-). See for
example:

James W. Hardin and Joseph M. Hilbe (2007) Generalized Linear Models
and Extensions, second edition. College Station: Stata Press.

or

 McCullagh, Peter; Nelder, John (1989). Generalized Linear Models,
Second Edition. Boca Raton: Chapman and Hall/CRC.

If you wish to explore those properties than I often find simulation
to be very illuminating, in addition to reading about those
techniques. Below is a quick set-up of how I would start with such a
simulation. For getting a reliable idea about the coverage I would use
more iterations: If the models works exactly as planned than we would
expect only 50 rejections of the true hypothesis out of a 1000
iterations, and it is the deviations from that number 50 that decides
how wrong or right the coverage is. With 20,000 iterations you would
expect a 1000 rejections of the true hypothesis, giving you much more
opportunity to detect deviations. Then the aim would be to explore how
well the model works and when it brakes down by trying out a set of
scenarios, e.g. different sample sizes.

*------------------ begin example -----------------------
program drop _all
program define sim, rclass	
	drop _all
	sysuse nlsw88
	gen y = exp(5 + .06*grade + 0.1*ttl_exp - ///
				0.002*ttl_exp*ttl_exp + rnormal(0, 1.041) )

	poisson y c.ttl_exp##c.ttl_exp grade, vce(robust)
	return scalar b_p  = _b[grade]
	return scalar se_p = _se[grade]
	glm y c.ttl_exp##c.ttl_exp grade, link(log) vce(robust)
	return scalar b_g  = _b[grade]
	return scalar se_g = _se[grade]
end

simulate b_p=r(b_p) se_p=r(se_p)   ///
         b_g=r(b_g) se_g=r(se_g) , ///
         reps(1000) : sim

simsum b*, se(se*) true(.06) mcse
*----------------- end example -----------------------------
(For more on examples I sent to the Statalist see:
http://www.maartenbuis.nl/example_faq )

This simulation requires Ian White's -simsum- package, which is
described in Ian White (2010) "simsum: Analyses of simulation studies
including Monte Carlo error" The Stata Journal, 10(3): 369-385. You
can install it by typing in Stata -ssc install simsum-.

Hope this helps,
Maarten

--------------------------
Maarten L. Buis
Institut fuer Soziologie
Universitaet Tuebingen
Wilhelmstrasse 36
72074 Tuebingen
Germany


http://www.maartenbuis.nl
--------------------------
*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/