[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

From |
"Nick Cox" <[email protected]> |

To |
<[email protected]> |

Subject |
st: RE: local macro in nlcom |

Date |
Tue, 18 Nov 2003 11:50:04 -0000 |

I don't understand what your problem is, but -nlcom- feeds on one or more formulae, not the result of a numerical evaluation. Here's an example of -nlcom- in action. I'm fitting a negative binomial distribution to a single response. This isn't necessarily your problem, but it may help in indicating some technique. These are real data: Hilborn and Mangel (1997, p.100) give data on incidental capture of albatrosses in the New Zealand subantarctic squid trawl fishery, 1990. (The birds get trapped accidentally in nets or trawl gear or cables.) Hilborn, R. and Mangel, M. 1997. The ecological detective: confronting models with data. Princeton University Press, Princeton, NJ. . input x freq 1. 0 807 2. 1 37 3. 2 27 4. 3 8 5. 4 4 6. 5 4 7. 6 1 8. 7 3 9. 8 1 10. 9 0 11. 10 0 12. 11 2 13. 12 1 14. 13 1 15. 14 0 16. 15 0 17. 16 0 18. 17 1 19. end Suppose I don't like the parameterisation used by -nbreg-. I can map to a new one, such as one commonly used in ecology. Here's the wrapper for -nlcom-: *! NJC 1.0.0 18 Nov 2003 program mynegbin version 8.0 syntax varname(numeric) [if] [in] [fweight] marksample touse di _n as txt "Fitting negative binomial distribution to `varlist'" qui nbreg `varlist' if `touse' [`weight' `exp'] nlcom (mean: exp(_b[_cons])) (k: exp(-_b[/lnalpha])) /// (p: exp(_b[_cons])/ (exp(_b[_cons]) + exp(-_b[/lnalpha]))) end Here are the results: . mynegbin x [fw=freq] Fitting negative binomial distribution to x mean: exp(_b[_cons]) k: exp(-_b[/lnalpha]) p: exp(_b[_cons])/ (exp(_b[_cons]) + exp(-_b[/lnalpha])) ---------------------------------------------------------------------- -------- x | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+-------------------------------------------------------- -------- mean | 2.777778 .2549476 10.90 0.000 2.27809 3.277466 k | 2.511678 .6363871 3.95 0.000 1.264382 3.758974 p | .5251538 .0672001 7.81 0.000 .393444 .6568636 ---------------------------------------------------------------------- -------- -nlcom- does the hard work for you. Incidentally, this little problem of mapping from one negative binomial parameterisation to another is also tackled in an ad hoc and Stata 6 way by -nbfit- on SSC. -nlcom- lets you do the job properly. Nick [email protected] Jun Xu > > Not sure if any solution or ideas about this problem. I try > to write some > routine for > getting the prediction and confidence interval for > predicted probabilities > in nbreg. However, > I got into some trouble. Note the following fourth line with: > > .local ai = "1/exp([lnalpha]_b[_cons])" // pay > attention to this > line > I think it is this line giving me trouble. I will get both > nonsensical > prediction and for sure incorrect standard errors. And if I have > > .local ai "1/exp([lnalpha]_b[_cons])" // pay > attention to this line > Again, I get exactly teh same as I did using the line above. > > If I take off the quotation mark;that is, if I have > > .local ai = 1/exp([lnalpha]_b[_cons]) // pay > attention to this line > then the prediction is CORRECT (I checked it), however, > because the local is > resolved before getting into nlcom, so the confidence > interval is incorrect > for sure. > > I got totally confused and don't know what to do. Also note > that `expxb' has > been defined to be equal to exp(xb)=mu and `i' refers to > the number of times > (y). Because I have tested it for many times and guess the > line that I just > talked about giving me trouble, especially how it is > translated in nlcom > affects the differences in the output, and I am pretty sure > it's not the > length of the expression list giving me trouble. Very > frustrated, and hope > someone could help out. In addition, I think the density > function I specify > for the nbreg is correct. > > ***************************************************************** > local alpha "exp([lnalpha]_b[_cons])" > *local alpha = "e(alpha)" > *local ai = 1/`alpha' > local ai = "1/exp([lnalpha]_b[_cons])" // > pay attention to > this line > local gai "exp(lngamma(`ai'))" > noi di in y `alpha' > noi di in y `ai' > noi di in y `expxb' > nlcom > ((exp(lngamma(`i'+`ai'))/(exp(lnfact(`i'))*exp(lngamma(`ai'))))* /// > ((`ai'/(`ai'+`expxb'))^`ai')* /// > ((`expxb'/(`ai'+`expxb'))^`i')) > ****************************************************************** > * * For searches and help try: * http://www.stata.com/support/faqs/res/findit.html * http://www.stata.com/support/statalist/faq * http://www.ats.ucla.edu/stat/stata/

**References**:**st: local macro in nlcom***From:*"Jun Xu" <[email protected]>

- Prev by Date:
**st: Re: the demise of cmdname and naming conventions** - Next by Date:
**Re: st: RE: Graphing a categorical variable: simple bar chart.** - Previous by thread:
**st: local macro in nlcom** - Next by thread:
**st: multiple line graph of long data** - Index(es):

© Copyright 1996–2024 StataCorp LLC | Terms of use | Privacy | Contact us | What's new | Site index |