On 31-Aug-07 13:06:42, Prof Brian Ripley wrote: > On Fri, 31 Aug 2007, Robin Hankin wrote: > >> Hi Kris >> lgamma() gives the log of the gamma function. > > Yes, but he used Igamma. According to ?pgamma, > > 'pgamma' is closely related to the incomplete gamma function. > As defined by Abramowitz and Stegun 6.5.1 > > P(a,x) = 1/Gamma(a) integral_0^x t^(a-1) exp(-t) dt > > P(a, x) is 'pgamma(x, a)'. Other authors (for example Karl > Pearson in his 1922 tables) omit the normalizing factor, > defining the incomplete gamma function as > 'pgamma(x, a) * gamma(a)'. > > and that seems to be what Igamma is following. GSL on the other > hand has the other tail, so > >> a <- 9 >> x <- 11.1 >> pgamma(x, a, lower=FALSE)*gamma(a) > [1] 9000.501 > >> You need gamma_inc() of the gsl package, a wrapper for the >> GSL library: >> >> > gamma_inc(9,11.1) >> [1] 9000.501 >> > > > As the above shows, you don't *need* this, but you do need the GSL > documentation to find out what R package gsl does. Why it differs from > the usual references is something for you to explain. Wikipedia > http://en.wikipedia.org/wiki/Incomplete_gamma_function > distinguishes them, as does MathWorld. > > I suggest you add a clarification to the gsl package as to what the > 'incomplete gamma function' means there.
We have been here before! -- though in connection with the Beta function in the first instance. See: See the thread starting on 13 Dec 2005 at http://finzi.psych.upenn.edu/R/Rhelp02a/archive/66670.html In particular I'll repeat my views on the distinction of terminology between "Incomplete Beta/Gamma Function" and "Beta/Gamma Distribution" (where "Function" refers to the incomplete *integral* and "Distribution" to the same divided by the complete integral i.e. by the "Beta/Gamma Function" which, in my view, should be defined as the complete integral). My reasons for preferring the terminology "Incomplete ... Function" for the incomplete integral *not* divided by the normalising constant (for both Beta and Gamma), and using "Distribution" for the incomplete integral divided by the constant (i.e. Pearson's "Ratio"), are several, but in summary: 1. The Beta and Gamma functions (not normalised) are fundamental mathematical functions in their own right; likewise their incomplete versions. 2. When needed in probability applications, then of course they need to be normalised; but then why not simply call them "distributions"? 3. (1) and (2) encapsulate in the terminology an essential distinction, and using (2) instead of (1) could lead to interesting inferences (e.g. that the complete Beta function is identically 1). I.e. the Beta function should not change its definition as x passes from 1 - epsilon to 1. And similarly for the Gamma. Granted there is non-uniformity of usage; but this does lead to confusion, which could be avoided by simply sticking to the distinction between "Incomplete ... Function" and "... Distribution". For more detail, see http://finzi.psych.upenn.edu/R/Rhelp02a/archive/66717.html (where, in particular, it is pointed out that both Karl Pearson and Abramowitz and Stegun are inconsistent, within the same publication, in their terminology, using "... Function" in one place to mean the integral, in another to mean the probability distribution. So it is unwise to appeal to either as the definitive reference, since the outcome will depend on where in the book you look it up). It looks as though the documentation for Igamma (ZipfR package) at http://finzi.psych.upenn.edu/R/library/zipfR/html/beta_gamma.html is admirably explicit as to how this (and related functions) are defined, so in this case there is no ambiguity. In the documenation for the Gamma functions in the gsl package, it is simply stated All functions [including gamma_inc()] as documented in the GSL reference manual section 7.19. There is no function named "gamma_inc" in the GSL reference manual. See: http://www.gnu.org/software/gsl/manual/html_node/Function-Index.html All functions are named like "gsl_sf_gamma_inc", so presumably this is what is intended; in which case it computes "the unnormalized incomplete Gamma Function \Gamma(a,x) = \int_x^\infty dt t^{a-1} \exp(-t) for a real and x >= 0." And again that is clear enough -- once you track it down! In many places in the R documentation (including the "?" pages) people have taken the trouble to spell out mathematical definitions (where these can be given in reasonable space). Especially in cases like the Incomplete Gamma and Beta functions, where there can be dispute over what is meant (see above), it is surely wise to spell it out! Best wishes to all, Ted. >> On 31 Aug 2007, at 00:29, [EMAIL PROTECTED] wrote: >> >>> Hello >>> >>> I am trying to evaluate an Incomplete gamma function >>> in R. Library Zipfr gives the Igamma function. From >>> Mathematica, I have: >>> >>> "Gamma[a, z] is the incomplete gamma function." >>> >>> In[16]: Gamma[9,11.1] >>> Out[16]: 9000.5 >>> >>> Trying the same in R, I get >>> >>>> Igamma(9,11.1) >>> [1] 31319.5 >>> OR >>>> Igamma(11.1,9) >>> [1] 1300998 >>> >>> I know I have to understand the theory and the math >>> behind it rather than just ask for help, but while I >>> am trying to do that (and only taking baby steps, I >>> must admit), I was hoping someone could help me out. >>> >>> Regard >>> >>> Kris. -------------------------------------------------------------------- E-Mail: (Ted Harding) <[EMAIL PROTECTED]> Fax-to-email: +44 (0)870 094 0861 Date: 01-Sep-07 Time: 15:49:57 ------------------------------ XFMail ------------------------------ ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.