Ken Ward's Mathematics Pages

Series: Stirling (for Gamma)

Page Contents

1. Stirling Formula Simple Proof
2. Statement of Further Approximations (Stirling Series)
3. Comparing Some Stirling Values with the Factorials

Stirling Formula Simple 'Proof'

We know from Euler's gamma function that the factorial of a number can be expressed as follows:

[1.01]

Where n≥0, and n is a real number.

Figure 1 below shows the graph of [1.01], when n=10. Its area is factorial ten, of course (because this is an exact formula, not an approximation),

An observation of the graph shows it is skewed to the left, and so is not symmetrical.

Let us call the curve we are integrating f(x), so:

[1.02]

By differentiating, equating to zero and solving, we obtain the value of n, for x. Differentiating again, gives a negative value at n (for positive n), indicating that when x=n, f(x) is maximum.

Let us transform the curve, f(x) by moving it n to the left. We do this so its maximum is at 0, and we are more able to compare f(x) with the Gaussian Integral.

In Figure 2, the cyan (light blue) curve is the Gaussian curve, and the fuchsia (pink) curve is the gamma(11). While we can see the two curves are similar, their absolute differences are quite great. For instance, at x=10, the gamma is about 23,000 at this point and the Gaussian is about 3,000. However, the maximum of both curves is nearly half a million, so the values at x=10 are relatively small.

(The Gaussian Integral appears here simply because I do not know how to integrate the gamma function, except for some values, and I wish eventually to use the Gaussian to find Stirling's formula.)

To transform the curve, let u=x+n and,  noting du=dx, we can substitute for x in [1.01], to obtain [1.03] below. Also, let us call the integral, I, for easy reference, so:

[1.03]

Because we have moved the curve, we have changed the limits on the integral, because when x=0, u=−n.

We can factorise (u+n)ⁿ by taking out the n giving n^n·(1+u/n)ⁿ, and separate the variable and constant parts of e (e^−(u+n)=e⁻u·e−n):

[1.04]

We can take the constants out of the integral, and we have:

[1.05]

We wish to simplify the integral by (eventually) exponentiating (1+u/n)ⁿ, so we take logarithms base e:

[1.06]

We can use Maclaurin's Series to expand the logarithm:

[1.07]

Hence, substituting this in [1.06], we obtain:

[1.08]

Now we can exponentiate [1.08], and substitute in [1.05], renaming the integral I₂. We have changed the limits (from ∞ to n), because our substitution requires |u|<n

[1.09]

We can now simplify the integral:

[1.10]

Noting that the integral:

[1.11]

is a Gaussian Integral with the indicated definite sum, but the limits of the Gaussian are infinity. That is, for I₂, [1.10], as n tends to infinity, I₂ tends I, which is n! So,

[1.12]

We therefore claim that as n tends to infinity, factorial n tends to the value below, which is Stirling's Approximation:

■  [1.13]

Statement of Further Approximations (Stirling Series)

It is possible to obtain greater accuracy in approximating factorials by using the Stirling Series, of which only three terms are given. The series can produce more accurate estimates up to a point, but after that it makes worse estimates. (This is stated without proof.)

[2.01]

Comparing Some Stirling Values with the Factorials

The following table compares some factorials with values calculated using Stirling's Approximation.

N	Factorial		%Error	Error Abs		%Error		%Error	Error Abs
1	1	0.92	7.79%	0.0779	1.00	0.10%	1.0025	-0.25%	-0.0025
2	2	1.92	4.05%	0.0810	2.00	0.05%	2.0007	-0.03%	-0.0007
3	6	5.84	2.73%	0.1638	6.00	0.03%	5.9989	0.02%	0.0011
4	24	23.51	2.06%	0.4938	24.00	0.02%	23.9960	0.02%	0.0040
5	120	118.02	1.65%	1.9808	119.99	0.01%	119.9862	0.01%	0.0138
6	720	710.08	1.38%	9.9218	719.94	0.01%	719.9404	0.01%	0.0596
7	5,040	4,980.40	1.18%	59.6042	5,039.69	0.01%	5,039.6863	0.01%	0.3137
8	40,320	39,902.40	1.04%	417.6045	40,318.05	0.00%	40,318.0454	0.00%	1.9546
9	362,880	359,536.87	0.92%	3,343.1272	362,865.92	0.00%	362,865.9180	0.00%	14.0820
10	3.6288E+06	3.5987E+06	0.83%	3.0104E+04	3.6287E+06	0.00%	3.6287E+06	0.00%	1.1525E+02
11	3.9917E+07	3.9616E+07	0.75%	3.0117E+05	3.9916E+07	0.00%	3.9916E+07	0.00%	1.0566E+03
12	4.7900E+08	4.7569E+08	0.69%	3.3141E+06	4.7899E+08	0.00%	4.7899E+08	0.00%	1.0728E+04
13	6.2270E+09	6.1872E+09	0.64%	3.9781E+07	6.2269E+09	0.00%	6.2269E+09	0.00%	1.1953E+05
14	8.7178E+10	8.6661E+10	0.59%	5.1729E+08	8.7177E+10	0.00%	8.7177E+10	0.00%	1.4502E+06
15	1.3077E+12	1.3004E+12	0.55%	7.2436E+09	1.3077E+12	0.00%	1.3077E+12	0.00%	1.9031E+07

From the graph we can notice that as the value of n increases, the relative (percentage ) error decreases, but the absolute error increases (at 15!, the error in the formula (without corrections) is more than 7 thousand million!)

We have argued that as n tends to infinity, the formula tends to n! However, even with small values of n, the approximations are quite close.