An approximation of the inverse of the gamma function is presented. The approximation is sufficiently accurate that it allows the inverse of the (discrete) factorial function to be expressed precisely.
This article is in response to my own previous inquiry (copied at the end), to which there were no responses, regarding the inverse of the gamma function.
Let k denote the positive zero of the digamma function, approximately 1.461632 . For x >= k, Gamma(x) is strictly increasing. Thus, restricting its domain accordingly, the inverse is a function. Perhaps this should be called the principal branch of the inverse of gamma. It is certainly the most interesting, and obtaining a good asymptotic approximation for it will be our main concern. We will eventually comment on the other branches; until then, it should be assumed that we are dealing with Gamma(x) for x >= k and its inverse.
Near the end of my recent sci.math article "23 kt. asymptotic expansion for gamma function", I presented a nice approximation for the gamma function:
Sqrt(2*pi)*((x-1/2)/e)^(x-1/2) - c
where c = Sqrt(2*pi)/e - Gamma(k), approximately 0.036534 .
This approximation of gamma is precisely invertible in terms of a well studied function; presumably, such cannot be said for the well know approximation Sqrt(2*pi/x)*(x/e)^x. The principal branch of the Lambert W function will be used. [Those unfamiliar with this function may consult
for details if necessary. Briefly, it is the inverse of the function x*e^x.] Letting L(x) = ln((x+c)/Sqrt(2*pi)), the inverse of my gamma approximation is
For large x, this gives a good approximation of the inverse of the gamma function. In the table below, I have chosen x = Gamma(N) = (N-1)! for various integer N purely for convenience, so that precise inverses are immediately obvious.
At the endpoint of the domain, we have AIG(Gamma(k)) = 1.50063, giving the approximation's worst |relative error| = 0.02668 . It should also be noted that the approximation's error itself, not just relative error, approaches 0 as x increases without bound. [Needless to say, such a strong statement cannot be made about the error in the approximation used here for the gamma function itself!]
A precise inverse for the factorial function N! for positive integer N can now be given by simply rounding to the nearest integer (since errors are sufficiently small -- less than 1/2 being all that was required):
InvFact(x) = Round(AIG(x)) - 1 for x = N!
As an example, InvFact(24) = Round(AIG(24)) - 1 = Round(4.995) - 1 = 5 - 1 = 4.
Asymptotic approximations for branches of the inverse of gamma other than the principal branch will now be considered briefly. For nonpositive integers n, let k_n denote the smallest zero of the digamma function greater than n. Then, for negative n, considering that Gamma(x) restricted to [ k_n, k_(n+1) ] is one-to-one, it has an inverse function, which we will call the n branch of the inverse of gamma; also, noting that k_0 = k, we separately define the 0 branch to be the previously discussed principal branch. (Apologies if my numbering of branches contravenes some well established system of which I am unaware!)
Nice _algebraic_ asymptotic approximations can be obtained fairly easily for all nonprincipal branches. (For all nonpositive n, the gamma function has a simple pole at n.) For example, for the -1 branch, the function 1/(g+x), where g is the Euler gamma constant (approximately 0.577216), works well when |x| is large.
The middle entry was included to show that such a simple approximation, not surprisingly, performs poorly near an endpoint of the domain. Algebraic approximations, involving square or cube roots, which perform far better near endpoints of the domain, as well as asymptotically, can also be derived fairly easily.
For all other branches, the simplest nice asymptotic approximations are rational functions having both numerator and denominator of degree 1. As examples: For the -2 branch, we have the approximation (2-g+x)/(g-1-x), and for the -3 branch, 4(g-2+2x)/(3-2g-4x).
Finally, as almost just a curiosity: Previously, for the Lambert W function, we used only its principal branch. However, if its other real branch is also used when obtaining the inverse of my gamma approximation, then AIG(x) = L(x) / W(L(x) / e) + 1/2 becomes bivalued close to the endpoint of its domain, now approximating part of the -1 branch of the inverse of gamma, as well as all of its principal branch. Although the approximation of that part of the -1 branch is poor, it is sufficiently accurate so that InvFact(x) = Round(AIG(x)) - 1 now gives the precise inverse relation of N! for integer N >= 0, with bivalued InvFact(1) being 0 or 1.
Comments are welcome.
David W. Cantrell
David W. Cantrell <DWCantrell@sigmaxi.org> wrote: > Anonymous <firstname.lastname@example.org> wrote: > > "Jeremy Price" <email@example.com> wrote in message > > news://Az6x7.103487$QK.firstname.lastname@example.org... > > > Yes, that was suppose to be the gamma function. > > > How do you find the inverse of it? > > > > Try a sketch of the gamma function, especially also on the -ve axis. > > Should give you some indication on part of the problem. However, even > > after restricting the domain, I don't think you can find a nice > > expression for the inverse, if that is what you are looking for. > > Related questions seem to occur fairly often in math newsgroups. > In my article at > , > I give, without proof, a possible inverse for the (discrete) factorial > function N! (for N >=1). > > I would be interested in knowing of results, including approximations, > concerning the inverse of gamma (with its domain restricted so that the > inverse is a function). Surely such an inverse has already been > investigated substantially. > > Regards, > David Cantrell