[R] Incomplete Gamma Functions and GammaDistribution Doc errata.

[Adrian Custer]

Hello all,

In the course of trying to implement the CDF of an
InverseGammaDistribution, I have run across the need for an igamma()
function. Several others have needed this function but the answers I
have found so far are not totally clear to me. I'm writing for three
reasons:
	1) to present a small error in the docs
	2) to clarify the approach we are expected to take
	3) to request,for the ease of use of myself and others, that the 	  R
developers implement igamma(), document the approach they 	  want users
to take or do both.

1) An error in the docs of the gamma distribution:
-------------------------------------------------

The GammaDistribution(param1,param2) is defined alternatively with two
versions of the second parameter which are inverses of one another.
'Rate' arises from one derivation such as that presented in Mathworld: 
http://mathworld.wolfram.com/GammaDistribution.html
(see just below equation (2) )
        "define $\theta \equiv 1 / \lambda$ to be the time between
        changes"
which makes \theta in equation (3) a rate parameter. In the R
documentation (i.e. ?pgamma) this same equation is presented but the
parameter \theta is now called s and is termed the scale. This goes
against general usage of mathworld, wikipaedia, the COLT library, and
others who call the "s" term a "rate". 

Please either change the equation or replace the name with "rate" and
the label with r or nu or theta.



2) Using existing functions to build the igamma(\alpha,\beta):
-------------------------------------------------------------

The incomplete gamma has been mentioned in two previous emails that I
could find:
http://www.r-project.org/nocvs/mail/r-help/2000/1006.html
and
http://www.r-project.org/nocvs/mail/r-help/2000/3618.html
both of which suggest that the upper incomplete gamma can be obtained 
using the pgamma() function. 

Is Prof. Ripley's message (the first above) saying that

	upperIncompletegamma(x,\rate) = pgamma(x,\rate) * gamma(\rate) ?

I have not yet been able to see why this would be true. I keep tripping
myself up on the differing notation of R, mathworld, wikipaedia and
Abramowitz and Stegun's "Handbook of mathematical functions".

Would someone please confirm that this is mathematically correct? 


3) Adding an igamma() function or documentation for how to generate it.
--------------------------------------------------------------------

Since I am at least the third person not clever enough to see the way to
derive the igamma directly, because the lack of an igamma() function has
cost me a chunk of my day, because using functions related to
distributions to derive pure mathematical functions is counter
intuitive, and because two implementations already exist, I hope the R
developers would be willing to incorporate an igamma() function directly
into R. 

The implementation could either use the code posted to the list:
http://www.r-project.org/nocvs/mail/r-help/1998/0295.html
or be a convinience function around the pgamma and gamma solution
suggested in the emails presented above. 

Alternatively, if this approach doesn't suit you, would someone please
add an explanation to the documentation for the gamma function
explaining how to generate an igamma(\alpha, \beta) using the current
functions?

Thanks to all,
adrian