Stirling’s approximation is often proved with heavy machinery, but there’s a short, elementary route starting from the integral \(n! = \int _0^\infty x^n e^{-x} \, dx\). With a single change of variables that turns \(\log (1+x) - x\) into \(-z^{2}/2\), the problem becomes Gaussian integration: the main mass contributes \(\sqrt {2\pi n}\), and the remainder is controlled by simple logarithmic inequalities. In a few lines this yields the sharp bounds \(\left (\frac {n}{e}\right )^n\sqrt {2\pi n} \le n! \le \left (\frac {n}{e}\right )^n\bigl (\sqrt {2\pi n}+1\bigr )\). The same viewpoint also exposes the higher-order corrections by integrating low-degree polynomials against the Gaussian, recovering the familiar \(1/(12n)\), \(1/(288n^2)\), \(\ldots \) terms; compare Robbins’ classical refinement of Stirling’s formula [3].
We will give a simple proof for the bound \[ (n/e)^n\sqrt {2\pi n} \le n! \le (n/e)^n(\sqrt {2\pi n}+1). \]
We use the integral representation of \(n!\): \begin {align*} n! &= \int _0^\infty x^n e^{-x} \\ &= n\int _0^\infty e^{n \log (nx) - nx} \, dx \\ &= n n^{n}\int _0^\infty e^{n (\log (x) - x)} \, dx \\ &= n (n/e)^{n}\int _{-1}^\infty e^{n (\log (1+x) - x)} \, dx. \end {align*}
Note that \(\log (1+x)-x\) is maximized at \(x=0\). Inspired by the series expansion \(\log (1+x)-x=-x^2/2 + O(x^3)\) we make a change of variables, \(\log (1+x)-x=-z^2/2\), or: \[ z = \mathrm {sign}(x)\,\sqrt {2\bigl (x-\log (1+x)\bigr )}. \] The differentials are simple: \[ d(-z^2/2)= -z\,dz = d(\log (1+x)-x) = \frac {dx}{1+x} - dx = -\frac {x}{1+x}dx, \] which after adjusting the integration limits yield: \[ n! = (n/e)^{n} \, n\int _{-\infty }^\infty e^{-nz^2/2} \, \bigl (z+z/x\bigr ) \, dz. \] We now need the bounds \begin {equation} 1+\tfrac {2}{3}z \;\le \; z+\frac {z}{x} \;\le \; \max \{1, 1+z\} \quad \text {for all }x\in \R . \label {eq:mid} \end {equation} Since \(\int _{-\infty }^\infty e^{-nz^2/2}dz=\sqrt {2\pi /n}\) and \(\int _{0}^\infty e^{-nz^2/2} z\, dz=1/n\), this immediately gives us \[ \sqrt {2\pi n} \le n\int _{-\infty }^\infty e^{-nz^2/2} \, (z+z/x) \, dz \le \sqrt {2\pi n} + 1. \] Which is what we want.
It remains to show \(\mathrm {(\ref {eq:mid})}\), which we will do using the following standard logarithmic inequalities [1]: \begin {align} \log (1+x) &\ge \frac {x(2+x)}{2+2x} \quad &&\text {when } x<0 \label {eq:logneg} \\ \log (1+x) &\ge x-\tfrac {1}{2}x^2 \quad &&\text {when } x>0 \label {eq:loglower} \\ \log (1+x) &\le \frac {x(6+x)}{6+4x} \quad &&\text {when } x>-1. \label {eq:logupper} \end {align}
We note \(z(1+1/x)\) is non-negative, since \(1+1/x<0\) implies \(\mathrm {sign}(x)=-1\). Thus it suffices to show \(z^2(1+1/x)^2\le 1\), which follows using \(\mathrm {(\ref {eq:logneg})}\): \[ z^2 = 2\bigl (x-\log (1+x)\bigr ) \le \frac {x^2}{1+x} \le \frac {x^2}{(1+x)^2} = (1+1/x)^{-2}, \] For \(x>0\) we can use \(\mathrm {(\ref {eq:loglower})}\) to bound \[ z = \sqrt {2\bigl (x-\log (1+x)\bigr )} \le \sqrt {2(x^2/2)} = x, \] so \(z/x\le 1\). This proves the upper bound of \(\mathrm {(\ref {eq:mid})}\). For the lower bound, it suffices to show \(z^2(1/3+1/x)^2 \ge 1.\) Using \(\mathrm {(\ref {eq:logupper})}\) we get \begin {align*} z^2\Bigl (\tfrac {1}{3}+\tfrac {1}{x}\Bigr )^2 &= 2\bigl (x - \log (1+x)\bigr )\Bigl (\tfrac {1}{3}+\tfrac {1}{x}\Bigr )^2 \\ &\ge 2\Bigl (x - \frac {x(6+x)}{6+4x}\Bigr )\Bigl (\tfrac {1}{3}+\tfrac {1}{x}\Bigr )^2 \\ &= \frac {x^2}{9 + 6 x} + 1 \quad \ge 1 \quad \text {since } x\ge -1. \end {align*}
The full Stirling approximation can similarly be derived from the series expansion: \[ z\Bigl (1+\frac {1}{x}\Bigr ) = 1 + \frac {2z}{3} + \frac {z^2}{12} - \frac {2z^3}{135} + \frac {z^4}{864} + \dots \] This form is nice as it allows easy integration with respect to the Gaussian pdf. E.g. \[ n\int _{-\infty }^{\infty }\frac {z^4}{864} \, e^{-n z^2/2} \, dz = \frac {\sqrt {2 \pi n}}{288n^{2}} \] matching the expected \[ n! \sim \sqrt {2\pi n}\left (\frac {n}{e}\right )^n \left (1 +\frac {1}{12n}+\frac {1}{288n^2} - \frac {139}{51840n^3} -\frac {571}{2488320n^4}+ \cdots \right ). \]
I think one can even give a uniform bound of \[ z\Bigl (1+\frac {1}{x}\Bigr ) \le 1 + \frac {2z}{3} + \frac {z^2}{c}, \] for some \(c\approx 8\), strengthening the upper bound even more. However, the logarithmic inequalities become quite tedious.
Finally, note that the change of variables gives \(x'(z) = \dfrac {dx}{dz} = z\bigl (1+1/x\bigr )\). From \(\log \bigl (z(1+1/x)\bigr ) \le 2z/3\) it follows that \(z\bigl (1+1/x\bigr ) \le e^{2z/3}\) for all \(z\). Plugging this into the integral representation \[ n! = \left (\tfrac {n}{e}\right )^n \; n \int _{-\infty }^{\infty } e^{-n z^2/2}\, \bigl (z+z/x\bigr )\, dz \] we obtain \[ n! \le \left (\tfrac {n}{e}\right )^n \; n \int _{-\infty }^{\infty } e^{-n z^2/2 + 2z/3}\, dz = \left (\tfrac {n}{e}\right )^n \sqrt {2\pi n}\, e^{2/(9n)}, \] where the last equality uses completion of the square. This yields the sharpened upper bound \(n! \le (n/e)^n\, \sqrt {2\pi n}\, e^{2/(9n)}\).
[1] Wikipedia contributors, “List of logarithmic identities — Inequalities,” Wikipedia, https://en.wikipedia.org/wiki/List_of_logarithmic_identities#Inequalities.
[2] T. D. Ahle, “Simple Stirling Upper/Lower Bounds,” MathOverflow answer, https://mathoverflow.net/a/484270/5429.
[3] H. Robbins, “A remark on Stirling’s formula,” The American Mathematical Monthly, 62(1):26–29, 1955. Available at https://heuklyd.github.io/papers/pdf/Robbins-1955.pdf.