[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]
st: Binomial regression
Constantine Daskalakis <C_Daskalakis@mail.jci.tju.edu>
st: Binomial regression
Thu, 02 Aug 2007 14:18:30 -0400
Good day to all.
I would like to share some observations on binomial regression using
Stata (and SAS).
Essentially, the problem is that we have a binary outcome Y (0/1) and want to model it as a function of covariates (X1, X2, etc) via the generalized linear model
p = b0 + b1*X1 + b2*X2 + ...
So far, I've done such projects in Stata and I've noted that binomial regression sometimes fails to converge and/or gives ML estimates of the coefficients that yield estimated probabilities that are outside the [0,1] range (typically, for some modest number of observations).
A colleague suggested that SAS might give better results and I was skeptical at first. But I had some time in my hands (it's summer after all!), so I took a few datasets and run them in Stata (-glm-) and SAS (Proc Genmod).
The Stata (8.2) implementation I used is the default ML (Newton-Raphson) algorithm
- glm y x1 x2 ..., fam(bin) link(i) search
- glm y x1 x2 ..., fam(bin) link(i) search fisher(#)
where # is some number of iterations with the Fisher scoring algorithm.
The SAS (9.1) implementation I used is
model y = x1 x2 ... / dist=bin link=id type3 wald itprint;
[I also played with specification of starting values, but I will leave this piece out, and confine my posting to results obtained using the default starting values in both packages.]
Here's what I've found:
Stata often gets bogged down ("backed up") after a few iterations and
does not converge.
Specifying Fisher scoring for some iterations in the beginning helps.
After Newton-Raphson takes over from Fisher scoring, it occasionally does converge. Most often, I have to use Fisher scoring throughout to get convergence. But see point #3 below.
SAS does seem to often converge (on the basis of parameter vector convergence), but also warns that the "relative Hessian convergence criterion" has not been achieved and that "convergence is questionable" (indicating that the likelihood has not really converged sufficiently).
(2) Likelihood of final model
The log-likelihood of the final Stata model is often somewhat better than that of the final SAS model. This might suggest that the Stata results are "better". However, see the drawback in point #4 below.
(3) Estimated coefficients and standard errors
Naturally, when Stata and SAS give different final models, their
estimated coefficients are different.
But beware using Fisher's scoring throughout to get convergence and a final model. Sometimes, this final model will have absurdly small standard errors (with p < 0.001 for all variables). If something like this happens, it might be useful to compute standard errors using the option "OPG":
- glm y x1 x2 ..., fam(bin) link(i) search fisher(#) opg
[There are special complications when there are covariate levels that have observed probability of 0 or 1 (ie, all observations are "0s" or "1s"), but I'll leave this issue aside.]
(4) Estimated probabilities
When Stata has convergence trouble (and sometimes when it does not), it
warns that some "parameter estimates produce inadmissible mean estimates
in one or more observations."
SAS gives no such warnings.
When I compute the predicted probabilities from both the SAS and the
Stata final models (which often are not the same, as I explained above),
SAS almost always has them in the [0,1] interval but Stata has quite a
few of them outside the interval.
I most recently ran 7 different models in a dataset (9 different outcomes, each with the same 4 predictor variables). My total N was 278.
For 2 outcomes, SAS and Stata gave identical results (with all predicted probabilities in the [0,1] range).
For the remaining 5 outcomes, SAS always gave a final model (although with the warning regarding the Hessian convergence).
Stata w/ NR did not converge after 50 iterations.
Stata w/ FS did converge in 4 of the 5 cases (after 6-20 iterations).
In 4 of the 5 cases, the final model of Stata had a somewhat better log-likelihood than the final model of SAS.
(i) Although SAS's estimated standard errors appeared sensible, Stata's (w/ FS) were clearly too small in 4 of the 5 cases.
(ii) The final SAS models never yielded a predicted probability less than 0 or greater than 1. In contrast, Stata's models yielded probabilities outside the [0,1] interval for 5-15% of the observations.
I've done this sort of comparison with a couple of other datasets with similar results. I know it is not proof positive (as opposed to a carefully constructed simulation), but I conclude that Stata's algorithm is not as robust as SAS's and would advise the use of SAS for binomial regression problems.
Perhaps I have missed something or perhaps other people have different experiences. So, please take this with the usual caution.
The documents accompanying this transmission may contain confidential
health or business information. This information is intended for the
use of the individual or entity named above. If you have received
this information in error, please notify the sender immediately and
arrange for the return or destruction of these documents.
Constantine Daskalakis, ScD
Thomas Jefferson University, Division of Biostatistics
1015 Chestnut St., Suite M100, Philadelphia, PA 19107
* For searches and help try: