Thai Journal of Mathematics

Volume 9 (2011) Number 1 : 21–28

www.math.science.cmu.ac.th/thaijournal
Online ISSN 1686-0209

More Accurate Approximations

for the Gamma Function

Gergő Nemes

Loránd Eötvös University, H-1117 Budapest,

Pázmány Péter sétány 1/C, Hungary
e-mail : [email protected]

Abstract : A series transformation idea inspired by a formula of R.W. Gosper

and some asymptotic expansions for the central binomial coefficients leads us to
new accurate approximations for the Gamma function.

Keywords : Asymptotic approximations; Asymptotic expansions; Gamma func-

tion; Laplace’s formula; Stirling’s formula.
2010 Mathematics Subject Classification : 33B15; 33F05; 41A60.

1 Introduction
The Gamma function plays an important role in several fields of mathematics
such as probability theory or combinatorics. One often has to evaluate the function
for large positive values. One way to aim this is to use asymptotic approximations.
It is well known that for large values of x the Gamma function has the asymptotic
series of the form [1, 2, 3]

√

x −x 1 1 139 571
Γ (x + 1) ∼ x e 2πx 1 + + − − + ··· .
12x 288x2 51840x3 2488320x4
(1.1)
Equation (1.1) is called Stirling’s formula however, Laplace was the first who
derived it by his approximation method for special integrals. Another famous

Copyright c 2011 by the Mathematical Association of Thailand. All

rights reserved.
22 T hai J. M ath. 9 (2011)/ G. Nemes

result is the Stirling series [1, 2, 3]

1 1 1 1 1
log Γ (x + 1) ∼ x + log x−x+ log (2π)+ − + −· · · . (1.2)
2 2 12x 360x3 1260x5

A main advantage of this latter series is that it has only odd powers of the variable.
For the past almost three hundred years several authors established fascinating new
asymptotic formulas to improve the accuracy of (1.1). For example, Karatsuba [4]
showed that a formula of Ramanujan can turn into an asymptotic expansion:
r
x −x √ 6 1 11
Γ (x + 1) ∼ x e π 8x3 + 4x2 + x + − + · · ·. (1.3)
30 240x
Mortici [5] proved in his more recent paper the following expansion similar to
Karatsuba’s:
r
n −n
√ 1 1 31 139
n! = Γ (n + 1) ∼ n e π 2n + + − − + · · ·. (1.4)
3 36n 3240n2 77760n3
For further developments in this topic, see, for example, [5, 6, 7, 8, 9, 10, 11, 12].
We develop some new variants of (1.1) in this paper and show that these
formulas are numerically more efficient than much of the early ones in many cases.
The first few values of the newly introduced coefficients and sequences can be
found in Appendix A.

2 New asymptotic expansions

The motivating examples are the following two asymptotic expansions [3, p.
12]:
22n

2n 1 1 5 21
∼ √ 1− + + − + · · · , (2.1)
n πn 8n 128n2 1024n3 32768n4
!
22n

2n 1 21
∼q 1− 2 + 4 − · · · . (2.2)
n π n + 14 64 n + 14 8192 n + 14

The first one is the standard asymptotic series of the central binomial coefficients.
If one expands them into a series in powers of 1/ (n + 1/4), the asymptotic series
contains only even powers. This remarkable result suggests that there might have
been an asymptotic expansion similar to (1.1) that is, it contains only even powers
of the shifted variable. The formula
s
x −x 1 1
Γ (x + 1) = x e 2π x + 1+O , (2.3)
6 x2

known as Gosper’s approximation [13] can be a good starting point. Our aim is
to elaborate the asymptotic series part in Gosper’s formula. It seems from (2.2)
More Accurate Approximations for the Gamma Function 23

that another series in terms of 1/ (x + 1/6) would be the right choice. It can be
shown that a series like that contains even and odd powers. If we insist to have
even powers only we are lead to the form
s !
x −x 1 g1 g2
Γ (x + 1) ∼ x e 2π x + g0 + 2 + 4 + ··· , (2.4)
6 (x + v1 ) (x + v2 )

where the sequences {gn }n≥0 and {vn }n≥1 has to be determined. One of our main
result is
Theorem 2.1. The Gamma function has an asymptotic series expansion of the
form  
s
1 X g n
Γ (x + 1) ∼ xx e−x 2π x + g0 +
2n
, (2.5)
6 (x + vn ) n≥1

as x → ∞, where the sequences {gn }n≥0 and {vn }n≥1 can be found from the
recurrence
n ⌊n/2⌋
−1/2 an−j −2j
gj vjn−2j , n ≥ 1.
X X
g0 = a0 = 1, j
= (2.6)
j=0
j 6 j=1
n − 2j

Here the an coefficients are those appearing in (1.1), i.e.,

Γ (x + 1) X an
√ ∼ . (2.7)
xx e−x 2πx n≥0 xn

Proof. As x → ∞ we have
−1/2 −1/2 X
Γ (x + 1) ex Γ (x + 1) ex

1 1 an
= √ 1+ ∼ 1+ .
xn
q
xx 2π x + 16
xx 2πx 6x 6x
n≥0

From the binomial formula we find

−1/2 X
1 −1/2 1
1+ ∼ ,
6x n 6n xn
n≥0

as x → ∞. Thus we obtain the asymptotic expansion

 
n
Γ (x + 1) ex

X X −1/2 a n−j  1 .
∼ (2.8)

j xn
q
xx 2π x + 1
j=0
j 6
6 n≥0

Suppose the expansion of the form

 
n
X X −1/2 a n−j  1 ∼ g0 +
X gn
2n . (2.9)

j 6 j xn
n≥0 j=0 n≥1 (x + vn )
24 T hai J. M ath. 9 (2011)/ G. Nemes

Now we expand the sum on the right-hand side in powers of 1/x:

vn −2n X gn X −2n vnl

X gn X gn
2n = 1+ =
(x + vn ) x2n x x2n l xl
n≥1 n≥1 n≥1 l≥0
 
X X −2n gn v l X ⌊n/2⌋
X −2j n−2j
= n
=  gj vj  1 .
l x2n+l j=1
n − 2j xn
n≥1 l≥0 n≥1

Clearly, from formula (2.9) we see that g0 = a0 = 1. What remains is to show

that the system

n ⌊n/2⌋
−1/2 an−j −2j
gj vjn−2j
X X
= (2.10)
j=0
j 6j j=1
n − 2j

has (unique) solutions {gn }n≥1 and {vn }n≥1 . Depending on the parity of n, we
have
2k k−1
X −1/2 a2k−j X −2j
gk = j
− gj vj2k−2j (2.11)
j=0
j 6 j=1
2k − 2j

and
2k+1
X k−1
−1/2 a2k−j+1 X −2j
− 2kgk vk = − gj vj2k−2j+1 , (2.12)
j=0
j 6j j=1
2k − 2j + 1

when k ≥ 1. First, expression (2.11) gives g1 = 1/144, and if we already found

the terms g1 , v1 , . . . , gk , then from formula (2.12) the value of vk follows. Now,
g1 , v1 , . . . , gk , vk determine gk+1 by (2.11). Thus, by induction, system (2.10)
defines the sequences {gn }n≥1 and {vn }n≥1 uniquely.

The first few values of the sequences {gn }n≥0 and {vn }n≥1 can be found in
Table 3. We have obtained an expansion in even powers however, the shift sequence
{vn }n≥1 has different terms whereas (2.2) has constant (= 1/4) shift in all terms.
Numerical evaluation of the first few vn (see Table 4) leads us to the

Conjecture 2.2. limn→∞ vn = 1/4.

This conjecture suggests a new asymptotic series to the Gamma function in

terms of 1/ (x + 1/4). Our second result is

Theorem 2.3. The Gamma function has an asymptotic series expansion of the
form
s
x −x 1 X Gn
Γ (x + 1) ∼ x e 2π x + n , (2.13)
6 x + 14
n≥0
More Accurate Approximations for the Gamma Function 25

as x → ∞, where the Gn coefficients are given by

n n−1
X
X −1/2 an−j −j Gj
Gn = − . (2.14)
j=0
j 6j j=0
n − j 4n−j

The an coefficients are from (1.1).

Proof. Similarly to the previous proof we expand our series in terms of 1/x:
X Gn X Gn X −n 1 1 X X −n Gn 1
n = =
x + 14 xn l 4l xl l 4l xn+l

n≥0 n≥0 l≥0 n≥0 l≥0
 
n
X X −j Gj  1
= .

j=0
n − j 4n−j  xn
n≥0

According to (2.8) and the uniqueness of asymptotic series the proof of (2.14) is
complete.

3 Numerical comparisons
We will compare in this paragraph the numerical performance of some asymp-
totic formulas to the Gamma function with our new formulas for large values. We
compare the following approximation formulas for Γ (x + 1).

x x √
1 1 1

2πx exp − + − ··· (Stirling) , (3.1)
e 12x 360x3 1260x5
x x √
1 1 139

2πx 1 + + 2
− − ··· (Laplace) , (3.2)
e 12x 288x 51840x3
s
x x √ 1 1 1

6
2πx 1+ + + − ··· (Ramanujan) , (3.3)
e 2x 8x2 240x3
s
x x √ 1 1 31

2πx 1+ + − − ··· (Mortici) , (3.4)
e 6x 72x2 6480x3
s !
x x 1 1 1
2π x + 1+ 2 − 3 − · · · (New) . (3.5)
e 6 144 x + 41 12960 x + 14

Table 1 displays the number of exact decimal digits (edd) of the formulas for
some values of x. Exact decimal digits are defined as follows:

approximation (x)
edd (x) = − log10 1 − . (3.6)
Γ (x + 1)
26 T hai J. M ath. 9 (2011)/ G. Nemes

In the table below the (i)-th entry (i = 1, 2, . . .) in a line starting with “name”
is the edd of the given approximation using the series up to the i-th order term.
The “−” sign indicates that the approximation is smaller and the “+” sign (not
displayed) indicates that the approximation is larger than the true value. Note
that in the case of Stirling’s formula the first n terms of the asymptotic series give
the 2nth order approximation.
Conclusion. It is seen that when we use odd order approximations, Laplace’s for-
mula is the most accurate. In the case of even orders Ramanujan’s approximation
is better than Stirling’s, Laplace’s and the one by Mortici, but our new formula
gives better approximations even than that of Ramanujan’s.

Formula x (1) (2) (3) (4) (5) (6) (7) (8)

Stirling 100 8.6 -13.1 17.2 -21.1
Laplace 100 -6.5 8.6 11.7 -13.1 16.2 17.2 -20.4 -21.1
Ramanujan 100 -5.7 -9.2 11.0 -13.3 -15.4 17.3 19.5 -21.1
Mortici 100 -6.2 8.6 11.4 -13.1 -15.9 17.2 -20.0 -21.1
New 100 -6.2 10.1 10.9 -14.9 -15.2 19.4 19.2 -23.0
Stirling 1000 11.6 -18.1 24.2 -30.1
Laplace 1000 -8.5 11.6 15.6 -18.1 22.2 24.2 -28.3 -30.1
Ramanujan 1000 -7.7 -12.2 15.0 -18.3 -21.4 24.3 27.5 -30.1
Mortici 1000 -8.2 11.6 15.4 -18.1 -21.9 24.2 -28.0 -30.1
New 1000 -8.2 13.1 14.9 -19.7 -21.2 26.9 27.2 -33.5
Stirling 10000 14.6 -23.1 31.2 -39.1
Laplace 10000 -10.5 14.6 19.6 -23.1 28.2 31.2 -36.3 -39.1
Ramanujan 10000 -9.7 -15.2 19.0 -23.3 -27.4 31.3 35.5 -39.1
Mortici 10000 -10.2 14.6 19.3 -23.1 -27.9 31.2 -36.0 -39.1
New 10000 -10.2 16.1 18.9 -24.7 -27.2 33.7 35.2 -42.1
Table 1: The number of exact decimal digits of the asymptotic series for
some values of x.
Expression (2.5) is a slightly different than the previous ones, thus we consider
it’s edds in a separate table. The notations are the same except the fact that this
series contains only even order terms.

Formula x (2) (4) (6) (8) (10)

Special 100 10.9 -15.2 19.2 -22.9 26.5
Special 1000 14.9 -21.2 27.2 -32.9 38.5
Special 10000 18.9 -27.2 35.2 -42.9 -50.5

Table 2: The number of exact decimal digits of the special asymptotic series
(2.5) for some values of x.
More Accurate Approximations for the Gamma Function 27

A Tables of coefficients

g0 1
1 23
g1 144 v1 90
3857 1792627
g2 − 3110400 v2 7289730
20932906335329 570984637359867601981
g3 34283052002304000 v3 2288928529497568067550

Table 3: The first few exact values of the sequences {gn }n≥0 and {vn }n≥1 .

g0 1.000000000000000
g1 0.006944444444444 v1 0.255555555555555
g2 −0.001240033436214 v2 0.245911302613402
g3 0.000610590513759 v3 0.249454987345193
g4 −0.000655407405149 v4 0.249839892410196
g5 0.001199164540953 v5 0.249958497082160

Table 4: The first few numerical values of the sequences {gn }n≥0 and
{vn }n≥1 .

53 360182239526821
G0 1 G5 − 2612736 G10 300361133850624000
5741173 104939254406053
G1 0 G6 9405849600 G11 210853515963138048000
1 37529 508096766056991140541
G2 144 G7 18811699200 G12 − 151814531493459394560000
1 710165119 70637580369737593
G3 − 12960 G8 − 1083553873920 G13 − 151814531493459394560000
257 3376971533 289375690552473442964467
G4 − 207360 G9 − 4022693756928000 G14 21861292535058152816640000

Table 5: The first few Gn coefficients.

Acknowledgements : I would like to thank C.J. Hegedűs and the anonymous

referee for their thorough, constructive and helpful comments and suggestions on
the manuscript.
28 T hai J. M ath. 9 (2011)/ G. Nemes

References
[1] M. Abramowitz, I.A. Stegun (eds.), Handbook of Mathematical Functions,
Dover Publications, 1965.
[2] E.T. Copson, Asymptotic Expansions, Cambridge University Press, 1965.
[3] Y.L. Luke, Mathematical Functions and their Approximations, Academic
Press, 1975.
[4] E.A. Karatsuba, On the asymptotic representation of the Euler gamma func-
tion by Ramanujan, J. Comp. Appl. Math. 135 (2001) 225–240.
[5] C. Mortici, Sharp inequalities related to Gosper’s formula, C. R. Acad. Sci.
Math. Paris 348 (2010) 137–140.
[6] C. Mortici, New approximation formulas for evaluating the ratio of gamma
functions, Math. Comput. Modelling 52 (2010) 425–433.
[7] C. Mortici, Asymptotic expansions of the generalized Stirling approximation,
Math. Comput. Modelling 52 (2010) 1867–1868.
[8] C. Mortici, The asymptotic series of the generalized Stirling formula, Comput.
Math. Appl. 60 (2010) 786–791.
[9] C. Mortici, New improvements of the Stirling formula, Appl. Math. Comput.
217 (2010) 699–704.
[10] C. Mortici, Ramanujan formula for the generalized Stirling approximation,
Appl. Math. Comput. 217 (2010) 2579–2585.
[11] G. Nemes, New asymptotic expansion for the Gamma function, Arch. Math.
(Basel) 95 (2010) 161–169.
[12] G. Nemes, On the coefficients of the asymptotic expansion of n!, J. Integer
Seqs. 13 (2010).
[13] R.W. Gosper, Decision procedure for indefinite hypergeometric summation,
Proc. Natl. Acad. Sci. 75 (1978) 40–42.

(Received 9 July 2010)

(Accepted 5 November 2010)

Thai J. Math. Online @ http://www.math.science.cmu.ac.th/thaijournal