How to calculate the probability Pr(y = n) after nbreg - math

[Alessandro Cucchetti]

After nbreg, what is the exact math during postestimation to obtain predicted probabilities when n=0, 1, 2.... k? This is not reported on the manual (or did I not found it?).

In other word, I want to replicate in a website the postestimation option pr(n) after nbreg
predict p_zero, p(0)
predict p_one, p(1)
predict p_two, p(2)
predict p_three, p(3)
predict p_four_or_more, p(4,.)

Thanks in advance to anyone who wants to help me.
Alessandro

[Alessandro Cucchetti]

Solved, in excel..

[John Mullahy]

Here's one way to do this in Stata, without resorting to Excel.

cap preserve
cap drop _all

set obs 1000
set seed 2345

gen x=.5+.5*runiform()
gen y=rnbinomial(1,x)

sum
loc ymax=r(max)

nbreg y x

forval j=0/`ymax' {
 predict p`j', pr(`j')
}

tab y
sum y p*

cap restore

Results:

. tab y

          y |      Freq.     Percent        Cum.
------------+-----------------------------------
          0 |        731       73.10       73.10
          1 |        172       17.20       90.30
          2 |         55        5.50       95.80
          3 |         31        3.10       98.90
          4 |          5        0.50       99.40
          5 |          2        0.20       99.60
          6 |          2        0.20       99.80
          8 |          1        0.10       99.90
         13 |          1        0.10      100.00
------------+-----------------------------------
      Total |      1,000      100.00

. sum y p*

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
           y |      1,000        .438    .9565935          0         13
          p0 |      1,000    .7320585    .1572274   .4075428   .9405516
          p1 |      1,000    .1718453    .0665825   .0559345   .2509211
          p2 |      1,000    .0557864    .0452106   .0033065   .1435028
          p3 |      1,000    .0219239    .0246726   .0001951   .0847285
-------------+---------------------------------------------------------
          p4 |      1,000    .0095707    .0132037   .0000115   .0499761
          p5 |      1,000    .0044463     .007116   6.77e-07     .02946
          p6 |      1,000    .0021514    .0038737   3.99e-08   .0173592
          p7 |      1,000    .0010715    .0021277   2.35e-09   .0102259
          p8 |      1,000    .0005453    .0011774   1.38e-10   .0060225
-------------+---------------------------------------------------------
          p9 |      1,000    .0002822    .0006555   8.13e-12   .0035463
         p10 |      1,000     .000148    .0003667   4.78e-13    .002088
         p11 |      1,000    .0000785     .000206   2.81e-14   .0012292
         p12 |      1,000     .000042    .0001161   1.65e-15   .0007236
         p13 |      1,000    .0000227    .0000657   9.73e-17   .0004259

(Added in edit) Of course you may wish to have predictions of outcome probabilities for outcome values outside the sample range. If so replace

forval j=0/`ymax'

with

forval j=0/`nn'

where nn is a local whose value is the maximum value for which you want a predicted probability.