Wrong confidence intervals using `margins` in Stata

Question

I am using Stata 13.1 to fit a logistic model and I am getting confidence intervals below 0 and above 1 when I predict probabilities using the margins command.

MRE:

sysuse auto, clear

* simple logistic regression
logit foreign mpg

* get predicted probabilities
margins, at(mpg=(5(5)40)) predict(pr)

* same result with expression
margins, at(mpg=(5(5)40)) exp(invlogit(predict(xb)))

The output is:

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         _at |
          1  |   .0271183   .0252542     1.07   0.283    -.0223792    .0766157
          2  |   .0583461   .0389888     1.50   0.135    -.0180704    .1347627
          3  |   .1210596   .0509373     2.38   0.017     .0212244    .2208948
          4  |   .2344013   .0547344     4.28   0.000      .127124    .3416787
          5  |   .4049667   .0743318     5.45   0.000      .259279    .5506543
          6  |   .6020462   .1162814     5.18   0.000     .3741387    .8299536
          7  |   .7707955   .1266899     6.08   0.000     .5224879    1.019103
          8  |   .8820117   .1004224     8.78   0.000     .6851874    1.078836

I'm trying to figure out what is Stata doing so I extracted the results with:

* save table as matrix
mat pr = r(table)
mat p = pr' //transpose matrix

Then I predicted the log odds (linear prediction), saved the results and transformed to probability (inverse logit)

* predict log odds
margins, at(mpg=(5(5)40)) exp(predict(xb))

* save table as matrix
mat tab = r(table)
mat t = tab' 

clear
svmat t

* transform logit to probability and display
gen prob   = invlogit(t1)
gen problo = invlogit(t5)
gen probhi = invlogit(t6)

format %9.3f prob*
list prob* in 1/8, noobs clean // correct confidence intervals?!

 prob   problo   probhi  
0.027    0.004    0.154  
0.058    0.015    0.199  
0.121    0.051    0.260  
0.234    0.144    0.358  
0.405    0.271    0.555  
0.602    0.369    0.797  
0.771    0.452    0.932  
0.882    0.530    0.980

Now if I take the results from the margins output in probability scale predict(pr) and (wrongly) use the SE from that output to produce CIs I get the same as Stata:

clear

* convert results to data
svmat p

* generate confidence intervals (the wrong way)
gen problo = p1 - (1.96*p2)
gen probhi = p1 + (1.96*p2)

format %9.3f p*
list p5 problo p6 probhi in 1/8, noobs clean

    p5   problo      p6   probhi  
-0.022   -0.022   0.077    0.077  
-0.018   -0.018   0.135    0.135  
 0.021    0.021   0.221    0.221  
 0.127    0.127   0.342    0.342  
 0.259    0.259   0.551    0.551  
 0.374    0.374   0.830    0.830  
 0.522    0.522   1.019    1.019  
 0.685    0.685   1.079    1.079

Sorry it took me so long but here is the question: Is there a way of getting the right CIs using the margins command directly? and does this problem happens with other GLMs?

Thanks

There's clearly a software component in this question, but I think we can keep this open since it goes deeper into the question of how to compute standard errors and confidence intervals for bounded range statistics (predicted probabilities). — StasK, Jul 31 '15 at 14:35

score 8 · Accepted Answer · answered Jul 31 '15 at 22:20

margins computes standard errors from nonlinear predictions using the delta-method, and as donlelek points out, it also uses a normal approximation for computing confidence intervals.

Consider a slight variation on donlelek's example.

sysuse auto, clear
logistic foreign mpg
margins, predict(pr) nopvalues

The result from margins is

Predictive margins                              Number of obs     =         74
Model VCE    : OIM

Expression   : Pr(foreign), predict()

--------------------------------------------------------------
             |            Delta-method
             |     Margin   Std. Err.     [95% Conf. Interval]
-------------+------------------------------------------------
       _cons |   .2972973   .0487662      .2017172    .3928773
--------------------------------------------------------------

Let's call this marginal prediction $m$.

Here we verify by hand the confidence limits produced by margins.

. display .2972973 - invnormal(.975)*.0487662
.2017173

. display .2972973 + invnormal(.975)*.0487662
.3928773

$m$ is the mean of the predicted probabilities. The predicted probability for observation $i$ is

$$ p_i = \frac{1}{1 + \exp(-z_i)} $$

where $z_i$ is the corresponding linear prediction

$$ z_i = b_0 + b_1*x_i $$

$x_i$ is the corresponding value of mpg, and $b$ is the vector of regression coefficients. Thus $m$ is

$$ m = \frac{1}{N}\sum_{i=1}^N p_i $$

$m$ is not a simple transformation of the marginal linear prediction. This is the general situation that margins was developed to handle. As such, the only method for computing a confidence interval for $m$ is the normal approximation using standard errors computed via the delta-method.

margins does not provide an option that will transform the confidence limits of a marginal linear prediction.

However, donlelek's question inspired me to write transform_margins, which will transform the point estimates and confidence limits from margins results with a user specified expression. Here is donlelek's example using transform_margins:

sysuse auto
logistic foreign mpg
margins, at(mpg=(5(5)40)) predict(xb)
transform_margins invlogit(@)

Here is the output from transform_margins

----------------------------------------------
             |         b         ll         ul
-------------+--------------------------------
         _at |
          1  |  .0271183   .0042517    .153952
          2  |  .0583461   .0151856   .1993467
          3  |  .1210596   .0511398   .2603462
          4  |  .2344013    .144129   .3575909
          5  |  .4049667   .2710298   .5547246
          6  |  .6020462   .3688264   .7966118
          7  |  .7707955   .4519781     .93203
          8  |  .8820117   .5300384   .9802168
----------------------------------------------

Type the following commands in Stata to install transform_margins.

net from http://www.stata.com/users/jpitblado
net describe transform_margins
net install transform_margins
help transform_margins

Thanks for your answer @Jeff Pitblado, just to clarify, there is no way to do this directly with margins, and the correct way of doing it would be to get the linear predictions and then transform to the desired scale with your transform_margins command. — donlelek, Aug 01 '15 at 23:45
margins does not have an option that provides the kind of confidence interval calculations you want. The normal approximation, using delta method standard errors, is the only option. In your example, I can see why you do not like this, so I provided transform_margins as a convenience. — Jeff Pitblado, Aug 03 '15 at 14:44
Thanks again @Jeff Pitblado, as an alternative; would it be ok to use the confidence intervals provided by margins and just: replace ll = 0 if ll < 0 and replace ll = 1 if ll > 1, and the same for ul? — donlelek, Aug 03 '15 at 14:51
I do not see why not, provided you document that the confidence intervals are based on a normal approximation and that you are cropping at the end-points. — Jeff Pitblado, Aug 03 '15 at 22:23

score 1 · Answer 2 · answered Jul 30 '15 at 21:01

1

Confidence intervals coming from predicted values in logistic regression models should be nonsymmetric and bounded between 0 and 1. The key is that the inverse link (the expit function) should be applied at the very end of calculation. Try the lincom command instead.

Margins uses marginal standardization which employs the delta method approximation to standard errors for predicted values on the probability scale.

answered Jul 30 '15 at 21:01

AdamO

62,637

Thanks for your answer @AdamO but my question is if there is a way of doing it right with margins. Also, this MRE is not my actual model and using lincom to do what I want would be tedious (I would if it's the only option) but it seems to me like a serious bug if Stata get the CIs wrong for all non-linear predictions. – donlelek Jul 31 '15 at 12:03

Wrong confidence intervals using `margins` in Stata

2 Answers2