Stirling's Formula - Approximation of n!

Posted on

I.Introduction: The Coin Tossing Problem

Whenever you have difficulty deciding something, for example, whether to attend the lectures of a course or not, the most common thing to do in this scenario is to toss a coin and see whether the heads or tails of the coin is facing up when they land. Assume that the coin you toss is uniform, meaning that the chance of the heads to land facing up is equal to the chance of the tails to land facing up. Then, the chance that you attend the lecture would be $\frac{1}{2}$. You could further explain that, if you always decide whether or not to attend a lecture by tossing a coin, you would probably end up attending about half of the lectures at the end of the semester (which is a relatively high attendance rate compared to somebody).

But then, what is the chance that you attend exactly half of the lectures in a semester. Mathematically put, if you have $2n$ lectures in total, what is the chance that you attend exactly $n$ lectures? There will be $2^{2n}$ possible outcomes in total if you toss the coin $2n$ times, and the possible outcomes that there are $n$ in $2n$ times when the coin faces up (you attend the lecture) will certainly be ${2n \choose n}$. Thus, the possibility that you attend exactly half of the lectures will be $$\begin{equation}P_{2n} = \frac{2n \choose n}{2^{2n}}=\frac{(2n)!}{2^{2n}n!n!}\label{ref1}\end{equation}$$

And when years and decades have passed, when $n\rightarrow\infty$ (assume that you would study in a university for the rest of your life), then, what would $P_{2n}$ in $\eqref{ref1}$ be? The problem immediately emerged is that how should we tackle the ugly $n!$ and $(2n)!$ ? A sensible way to deal with this is to approximate $n!$ with other terms, and this is when the powerful Stirling’s formula shall be introduced.

II.The Proof: Stirling’s Formula

Before getting our hands dirty into mathematical statements and equations, let us first take a glimpse and see how the formula looks like $$\begin{equation}\lim_{n\to\infty}\frac{n!}{\sqrt{2\pi}\cdot n^{n+\frac{1}{2}}\cdot e^{-n}}=1\label{ref2}\end{equation}$$

Quite intimidating, isn’t it? All the famous irrational numbers such as $\sqrt{2}$, $\pi$, and $e$ play a part in it. But if you’ve studied hard enough in your calculus course (at least attend more than half of the calculus lectures), the proof of this formula will be just a piece of cake for you. Don’t believe me? Let me give a proof of that!

It is quite obvious from the formula that if we take the natural log on both sides, it will be easier to handle, since many of the terms include exponents and there is an $e$ being the base number in one of the terms as well. So, taking $ln$ on both sides, and we’ll have $$\lim_{n\to\infty}ln\frac{n!}{\sqrt{2\pi}\cdot n^{n+\frac{1}{2}}\cdot e^{-n}}=ln\, 1$$

After arranging it properly using techniques of logarithm, we will further get $$\begin{equation}\lim_{n\to\infty}(ln\, n!-\frac{1}{2}ln\, 2-\frac{1}{2}ln\, \pi- n\, ln\, n- \frac{1}{2}ln\, n+n)=0\label{ref3}\end{equation}$$

This should seem less intimidating. Our mission now becomes approximating the term $ln\, n!$ using $-\frac{1}{2}ln\, 2-\frac{1}{2}ln\,\pi-n\, ln\, n- \frac{1}{2}ln\, n + n$. We’ll start by showing that $$\begin{equation}ln\, n!=\sum_{k=1}^{n-1}\int_{k-\frac{1}{2}}^{k+\frac{1}{2}}(ln\, k)dt + \int_{n-\frac{1}{2}}^{n}(ln\, n)dt+\frac{1}{2}ln\, n\label{ref4}\end{equation}$$

By calculating the integrations in the RHS, we’ll get

$$\begin{align} RHS & = \int_{1-\frac{1}{2}}^{1+\frac{1}{2}}(ln\, 1)dt+\int_{2-\frac{1}{2}}^{2+\frac{1}{2}}(ln\, 2)dt+ \dotsb +\int_{n-1-\frac{1}{2}}^{n-1+\frac{1}{2}}(ln(n-1))dt+\int_{n-\frac{1}{2}}^{n}(ln\, n)dt+\frac{1}{2}ln\, n \newline & = (ln\, 1)+(ln\, 2)+\dotsb +[ln(n-1)]+\int_{n-\frac{1}{2}}^{n}(ln\, n)dt+\frac{1}{2}ln\, n \newline & = ln(n-1)!+\frac{1}{2}ln\, n+\frac{1}{2}ln\, n \newline & = ln\, n!\end{align}$$

We have the identity that

$$I(x)=\int_{0}^{x}(ln\, t)dt=x\, ln\, x - x$$

so it is natural to approximate the sum in $\eqref{ref4}$ by

$$\sum_{k=1}^{n-1}\int_{k-\frac{1}{2}}^{k+\frac{1}{2}}(ln\, t)dt+\int_{n-\frac{1}{2}}^{n}(ln\, t)+\frac{1}{2}ln\, n=\int_{\frac{1}{2}}^{n}(ln\, t)dt+\frac{1}{2}ln\, n = I(n)-I(\frac{1}{2})+\frac{1}{2}ln\, n$$

Hence, what remains is to bound the error of the approximation. To do so, let us define that

$$\begin{align}a_k & = \frac{1}{2}ln\, k-\int_{k-\frac{1}{2}}^{k}(ln\, t)dt \newline b_k & = \int_{k}^{k+\frac{1}{2}}(ln\, t)dt-\frac{1}{2}ln\, k\end{align}$$

where $k = 1, 2, \dotsb$.

As $a_k-b_k=ln\, k-\int_{k-\frac{1}{2}}^{k+\frac{1}{2}}(ln\, t)dt$,

$$\begin{equation}\begin{split}\sum_{k=1}^{n-1}(a_k-b_k)+a_n & = \sum_{k=1}^{n-1}ln\, k-\sum_{k=1}^{n-1}\int_{k-\frac{1}{2}}^{k+\frac{1}{2}}ln\, t\, dt+\frac{1}{2}ln\, n-\int_{n-\frac{1}{2}}^{n}ln\, t\, dt \newline & = ln\, n!-\frac{1}{2}ln\, n-\int_{\frac{1}{2}}^{n}ln\, t\, dt \newline & = ln\, n!-\frac{1}{2}ln\, n-I(n)+I(\frac{1}{2})\end{split}\label{ref5}\end{equation}$$

The result of $\eqref{ref5}$ implies that the alternating sum converges to a finite number (i.e. $S=\sum_{k=1}^{\infty}(a_k-b_k) < \infty$), which is a helpful statement that we’ll make use of further in the proof. But before moving onto that, we’ll first develop some other useful tools required to complete the proof, that is, for $k=1,\, 2,\,\dotsb$,



This could also be proven trivially by integration,

$$\begin{equation}\begin{split}a_k & = \frac{1}{2}ln\, k-\int_{k-\frac{1}{2}}^{k}(ln\, t) dt\qquad(by\,definition) \newline & = \frac{1}{2}[k\, ln\, k-k-(k -\frac{1}{2})ln(k-\frac{1}{2})+k-\frac{1}{2}] \newline & = \frac{1}{2}+\frac{1}{2}ln\, k- k\, ln\, k+(k-\frac{1}{2})ln(k-\frac{1}{2}) \newline & = \frac{1}{2}-(k-\frac{1}{2})ln\, k+(k - \frac{1}{2})ln(k-\frac{1}{2}) \newline & = \frac{1}{2}-(k-\frac{1}{2})ln \frac{k}{k-\frac{1}{2}}\end{split}\label{ref8}\end{equation}$$

Looking back at $\eqref{ref6}$,

$$\begin{split}RHS & = \frac{1}{2}ln(\frac{1}{1-\frac{1}{2k}})+k\, ln(k-\frac{1}{2})+\frac{1}{2}-k\, ln\, k \newline & = \frac{1}{2}-k\, ln \frac{k}{k-\frac{1}{2}}+\frac{1}{2}ln \frac{1}{1-\frac{1}{2k}} \newline & = \frac{1}{2}-k\, ln \frac{k}{k-\frac{1}{2}}+\frac{1}{2}ln\frac{k}{k-\frac{1}{2}} \newline & = \frac{1}{2}-(k-\frac{1}{2})ln \frac{k}{k-\frac{1}{2}} \newline & = a_k\qquad by\,\eqref{ref8}\end{split}$$

Similarly, for $b_k$,

$$\begin{equation}\begin{split}b_k & = \int_{k}^{k+\frac{1}{2}}(ln\, t)dt-\frac{1}{2}ln\, k\quad (by\, definition) \newline & = (k+\frac{1}{2})ln(k+\frac{1}{2})-(k+\frac{1}{2})-k\, ln\, k+k-\frac{1}{2}ln\, k\newline & = -\frac{1}{2}+(k+\frac{1}{2})ln(k+\frac{1}{2})-(k+\frac{1}{2})ln\, k \newline & = -\frac{1}{2}+(k+\frac{1}{2})ln(\frac{k+\frac{1}{2}}{k})\end{split}\label{ref9}\end{equation}$$

Looking back at $\eqref{ref7}$, $$\begin{split}RHS & = -\frac{1}{2}+\frac{1}{2}ln(\frac{1}{2k}+1)+k\, ln(\frac{1}{2}+k)-k\, ln\, k\newline & = -\frac{1}{2}+\frac{1}{2}ln(\frac{1}{2k}+1)+k\, ln (\frac{k+\frac{1}{2}}{k}) \newline & = -\frac{1}{2}+\frac{1}{2}ln(\frac{k+\frac{1}{2}}{k})+k\, ln(\frac{k+\frac{1}{2}}{k})\newline & = -\frac{1}{2}+(k+\frac{1}{2})ln(\frac{k+\frac{1}{2}}{k})\qquad by\,\eqref{ref9}\end{split}$$

We’ve already discussed that $S=\sum_{k=1}^{\infty}(a_k-b_k)<\infty$ after proving $\eqref{ref5}$, now that we have proven $\eqref{ref6}$ and $\eqref{ref7}$, we’re finally able to figure out what $S$ is. In order to prove the Stirling’s formula, we should write $S$ as $$\begin{equation}S=I(\frac{1}{2})-\int_{0}^{\frac{1}{2}}ln\, \phi (t)dt\label{ref10}\end{equation}$$

where $$\phi (t)=t\prod_{k=1}^{\infty}(1-\frac{t^2}{k^2})=\frac{\sin\,\pi t}{\pi}$$

Let’s prove it by using the identities of $\eqref{ref6}$ and $\eqref{ref7}$,

$$\begin{align}S & = \sum_{k=1}^{\infty}(a_k-b_k)\qquad (by\, definition) \newline & = \sum_{k=1}^{\infty} ( \int_{0}^{\frac{1}{2}}ln\frac{1}{1-\frac{t}{k}}dt-\int_{0}^{\frac{1}{2}}ln(1+\frac{t}{k})dt ) \qquad (by\,\eqref{ref6},\,\eqref{ref7}) \newline & = \sum_{k=1}^{\infty}\int_{0}^{\frac{1}{2}}(ln\frac{k}{k-t}-ln\frac{k+t}{k})dt \newline & = \sum_{k=1}^{\infty}\int_{0}^{\frac{1}{2}}ln\frac{k\cdot k}{(k-t)(k+t)}dt \newline & = \sum_{k=1}^{\infty}\int_{0}^{\frac{1}{2}}ln\frac{k^2}{k^2-t^2}dt \newline & = -\sum_{k=1}^{\infty}\int_{0}^{\frac{1}{2}}ln\frac{k^2-t^2}{k^2}dt \newline & = -\sum_{k=1}^{\infty}\int_{0}^{\frac{1}{2}}ln(1-\frac{t^2}{k^2})dt \newline & = -\int_{0}^{\frac{1}{2}} (ln(1-\frac{t^2}{1^2})+ln(1-\frac{t^2}{2^2})+\dotsb +ln(1-\frac{t^2}{n^2}))dt\quad (n\to\infty) \newline & = -\int_{0}^{\frac{1}{2}}ln(\prod_{k=1}^{\infty}(1-\frac{t^2}{k^2}))dt \newline & = \int_{0}^{\frac{1}{2}}ln\, t\, dt-\int_{0}^{\frac{1}{2}}ln\, t\, dt-\int_{0}^{\frac{1}{2}}ln(\prod_{k=1}^{\infty}(1-\frac{t^2}{k^2}))dt \newline & = \int_{0}^{\frac{1}{2}}ln\, t\, dt-\int_{0}^{\frac{1}{2}}ln\, t(\prod_{k=1}^{\infty}(1-\frac{t^2}{k^2}))dt \newline & = I(\frac{1}{2})-\int_{0}^{\frac{1}{2}}ln\,\phi (t)\, dt\end{align}$$

By the fact that $\phi (t)=t\prod_{k=1}^{\infty}(1-\frac{t^2}{k^2})=\frac{\sin\,\pi t}{\pi}$, we could arrange $\eqref{ref10}$ and rewrite it as $$\begin{equation}S-I(\frac{1}{2})=\frac{1}{2}ln\,\pi-\int_{0}^{\frac{1}{2}}ln(\sin\,\pi t)dt\label{ref11}\end{equation}$$

The term $-\int_{0}^{\frac{1}{2}}ln(\sin\,\pi t)dt$ seems annoying, but fortunately, it could be transformed into a simpler term that we’ve seen before.

Since $$\begin{split}\int_{0}^{\frac{1}{2}}ln(\sin\,\pi t)dt & = \int_{0}^{\frac{1}{2}}ln(\cos(\pi(\frac{1}{2}-t)))dt \newline & = \int_{0}^{\frac{1}{2}}ln(\cos\,\pi u)du\qquad (by\, letting\, u=\frac{1}{2}-t) \newline & = \int_{0}^{\frac{1}{2}}ln(\cos\,\pi t)dt\end{split}$$


$$\begin{split}\int_{0}^{\frac{1}{2}}ln(\sin\, 2\pi t)dt & = \frac{1}{2}\int_{0}^{1}ln(\sin\,\pi u)du\qquad (by\, letting\, u=2t) \newline & = \int_{0}^{\frac{1}{2}}ln(\sin\,\pi u)du\qquad (by\, symmetry)\newline & = \int_{0}^{\frac{1}{2}}ln(\sin\,\pi t)dt\end{split}$$

we could get

$$\begin{split}\int_{0}^{\frac{1}{2}}ln(\sin\,\pi t)dt & = \int_{0}^{\frac{1}{2}}ln(\sin(\pi t))dt+\int_{0}^{\frac{1}{2}}ln(\cos(\pi t))dt-\int_{0}^{\frac{1}{2}}ln(\sin(2\pi t))dt \newline & = \int_{0}^{\frac{1}{2}}ln\frac{\sin (\pi t)\,\cos (\pi t)}{\sin (2\pi t)}dt \newline & = \int_{0}^{\frac{1}{2}}ln\frac{1}{2}dt\qquad (by\, double\, angle\, formula)\newline & = -\frac{1}{2}ln\, 2\end{split}$$

Thus, the equation $\eqref{ref11}$ could be beautifully transformed into $$\begin{equation}S-I(\frac{1}{2})=\frac{1}{2}ln\,\pi+\frac{1}{2}ln\, 2\label{ref12}\end{equation}$$

Now we’ve finally got every puzzle pieces, let’s put them together to prove Stirling’s formula. The basic idea is to merge the two $S-I(\frac{1}{2})$ we know together.

$$\begin{split}S-I(\frac{1}{2}) & = \lim_{n\to\infty}\sum_{k=1}^{n}(a_k-b_k)-I(\frac{1}{2})\qquad (by\, definition) \newline & = \lim_{n\to +\infty}\sum_{k=1}^{n-1}(a_k-b_k)+a_n-I(\frac{1}{2}) \newline & \quad(\because\lim_{n\to +\infty}b_n=\lim_{n\to +\infty}\int_{0}^{\frac{1}{2}}ln(1+\frac{t}{n})dt=0) \newline & = \lim_{n\to +\infty}(ln\, n!-\frac{1}{2}ln\, n-I(n))\qquad (by\,\eqref{ref5})\end{split}$$

We could also write $S-I(\frac{1}{2})$ as $$S-(\frac{1}{2})=\frac{1}{2}ln\,\pi+\frac{1}{2}ln\, 2\qquad\eqref{ref12}$$

So, $$\lim_{n\to +\infty}[ln\, n!-\frac{1}{2}ln\, n-I(n)]=\frac{1}{2}ln\,\pi+\frac{1}{2}ln\, 2$$

After arranging, we have

$$\lim_{n\to +\infty}[ln\, n!-\frac{1}{2}ln\, n-I(n)]-\frac{1}{2}ln\,\pi-\frac{1}{2}\ln\, 2=0$$

$$\lim_{n\to\infty}(ln\, n!-\frac{1}{2}ln\, 2-\frac{1}{2}ln\,\pi-n\, ln\, n-\frac{1}{2}ln\, n+n)=0$$

which is the same as $\eqref{ref3}$ and the original $\eqref{ref2}$, $$\lim_{n\to\infty}\frac{n!}{\sqrt{2\pi}\cdot n^{n+\frac{1}{2}}\cdot e^{-n}}=1$$

Congratulations! We’ve completed the proof, now it’s time for us to put it back into our original coin tossing problem, that is, $$P_{2n} = \frac{2n \choose n}{2^{2n}}=\frac{(2n)!}{2^{2n}n!n!}\qquad\qquad\eqref{ref1}$$

Substitute the $(2n)!$ and $n!$s with what we acquired in the Stirling’s formula, $$P_{2n}=\frac{\sqrt{2\pi}\cdot (2n)^{2n+\frac{1}{2}}\cdot e^{-2n}}{2^{2n}\cdot (\sqrt{2\pi}\cdot n^{n+\frac{1}{2}}\cdot e^{-n})^{2}}$$

Arrange it, and we’ll have

$$P_{2n}=\frac{1}{\sqrt{\pi n}}$$

When $n\to\infty$,

$$P_{2n}=\lim_{n\to +\infty}\frac{1}{\sqrt{\pi n}}=0$$

Surprisingly, it equals to zero, meaning that it is very unlikely that you would attend exactly half of the lectures after many days has passed.

III. A More Quantitative formula: Bounding the error

We have now seen the proof and a simple application of Stirling’s formula. However, since it is a method by approximating the value of $n!$, we had better know the bounds of the error in order to fully understand the formula. he more quantitative error bounding formula is simpy $$1<\frac{n!}{\sqrt{2\pi\cdot n^{n+\frac{1}{2}}\cdot e^{-n}}}<e^{\frac{1}{8n}}$$ for any $n\geq 1$.

The proof of this equation requires some of the tools we developed earlier. From $\eqref{ref6}$ and $\eqref{ref7}$, we know that $a_k>b_k>a_{k+1}>0$ because for $t\in (0, \frac{1}{2} ] $,



We have $b_k>0$ for all $k$.

Also, $$a_k-b_k=\int_{0}^{\frac{1}{2}}ln\frac{k^2}{k^2-t^2}dt\qquad (from\,\eqref{ref6},\,\eqref{ref7})$$

for $t\in (0, \frac{1}{2}]$,



We have $a_k-b_k>0$ for all $k$.

And that,

$$b_k-a_{k+1}=\int_{0}^{\frac{1}{2}}ln\frac{(k^2+k)-(t^2-t)}{k^2+k}dt\qquad (from\,\eqref{ref6},\,\eqref{ref7})$$

for $t\in (0,\frac{1}{2}]$,



We have $b_k-a_{k+1}>0$ for all $k$.

Combine the above three inequalities, we have


Let us further define the approximation error by

$$E_n=ln\, n!-\frac{1}{2}ln\, n-I(n)-[S-I(\frac{1}{2})]=b_n-a_{n+1}+b_{n+2}-\dotsb$$

By $\eqref{ref13}$, we have

$$E_n=(b_n-a_{n+1})+(b_{n+1}-a_{n+2})+\dotsb >0$$

$$b_n-E_n=(a_{n+1}-b_{n+1})+(a_{n+2}-b_{n+2})+\dotsb >0$$

From this result, we could show that $$\begin{equation}0<E_n<b_n\leq\frac{1}{8n}\label{ref14}\end{equation}$$


$$\begin{split}b_n & = \int_{0}^{\frac{1}{2}}ln(1+\frac{t}{n})dt\qquad (from\,\eqref{ref7}) \newline \leq & \int_{0}^{\frac{1}{2}}\frac{t}{n}dt \newline \qquad &(\because\, ln(1+x)\leq x\, for\, x \geq 0) \ =&\frac{1}{8n}\end{split}$$

Transform $\eqref{ref14}$ into $$e^0 < e^{E_n} < e^{\frac{1}{8n}}$$

Plug $E_n=ln\, n!-\frac{1}{2}ln\, n-I(n)-[S-I(\frac{1}{2})]$ and $e^0 = 1$ into it, we’ll get:

$$1<e^{(ln\, n!-\frac{1}{2}ln\, n-I(n)-[S-I(\frac{1}{2})])}<e^{\frac{1}{8n}}$$

After arranging, we have $$1<\frac{n!}{\sqrt{2\pi\cdot n^{n+\frac{1}{2}}\cdot e^{-n}}}<e^{\frac{1}{8n}}$$

which is the desired error bound.

IV. Conclusion: Further Applications

After understanding the proof of Stirling’s formula, it’s time for us to discover the applications of it.

Generally, whenever you have to approximate $n!$, Stirling’s formula would very likely be the tool you need to use. For example, when calculating probabilities, we often run into the problematic term $n!$, just like what we have encountered in the first example.

Another place where Stirling’s formula is useful is in the field of computer science. When you write a function to calculate the value of $n!$, the most intuitive way is to multiple all the numbers from $2$ to $n$ to get the final result, just as the following implementation.

{{ highlight python linenos=“table”}} factorial(n): if (n=0) return 1 else return n*factorial(n-1) {{ / highlight}}

This algorithm is linear, meaning that it would run in $O(n)$ time. However, if you write the function with another algorithm, using Stirling’s formula,

{{ highlight python linenos=“table” }} factorial(n): return sqrt(PI*2)*pow(n / E, n) {{ / highlight }}

it would be less time-consuming, since the calculation of $n$ to the $n$-th power (aka. the $(\frac{n}{e})^{n}$ is only $O(log\, n)$.

Through this simple yet elegant equation, much of the complexities in algorithms could be reduced and let us understand further in the fields of mathematics and engineering.