史特灵公式

史特灵公式（英语：Stirling's formula）是一条用来取n阶乘近似值的数学公式。一般来说，当n很大的时候，n阶乘的计算量十分大，所以史特灵公式十分好用，而且，即使在n很小的时候，史特灵公式的取值已经十分准确。这个公式以詹姆斯·史特灵（英语：James Stirling (mathematician)）的名字命名，虽然亚伯拉罕·棣美弗早于史特灵提出了一个类似的公式，但结果较不精确。^[1][2][3]

当n增加时，(ln n!)与o (n ln n − n)之比趋于1

阶乘与史特灵公式

${\displaystyle n!}$（蓝色）、${\displaystyle {\sqrt {2\pi n}}\left({\frac {n}{e}}\right)^{n}}$（橘色），数字越大${\displaystyle {\sqrt {2\pi n}}\left({\frac {n}{e}}\right)^{n},}$会越趋近${\displaystyle n!}$。但${\displaystyle {\sqrt {2\pi n}}\left({\frac {n}{e}}\right)^{n}}$在负值则会因为出现虚数而无法使用

史特灵公式为：

${\displaystyle n!\approx {\sqrt {2\pi n}}\,\left({\frac {n}{e}}\right)^{n}.}$

这就是说，对于足够大的整数n，这两个数互为近似值。更加精确地：

${\displaystyle \lim _{n\rightarrow \infty }{\frac {n!}{{\sqrt {2\pi n}}\,\left({\frac {n}{e}}\right)^{n}}}=1}$

或

${\displaystyle \lim _{n\rightarrow \infty }{\frac {e^{n}\,n!}{n^{n}{\sqrt {n}}}}={\sqrt {2\pi }}.}$

历史

这个公式是亚伯拉罕·棣莫弗首先发现的，形式为：

${\displaystyle \displaystyle n!\approx cn^{n+{\frac {1}{2}}}e^{-n}}$，其中c为常数。

斯特林（英语：James_Stirling_(mathematician)）证明了公式中${\displaystyle c={\sqrt {2\pi }}}$，约为2.506628274631。更加精确的形式是雅克·比内发现的。

推导

方法一

这个公式，以及误差的估计，可以推导如下。首先不直接估计n!，而是考虑它的自然对数：

${\displaystyle \ln(n!)=\ln 1+\ln 2+\cdots +\ln n.}$

即:

${\displaystyle \ln(n!)-{\frac {\ln n}{2}}=\ln 1+\ln 2+\cdots +\ln n-{\frac {\ln n}{2}}.}$

这个方程的右面是积分${\displaystyle \int _{1}^{n}\ln(x)\,dx=n\ln n-n+1}$的近似值（利用梯形法则），而它的误差由欧拉-麦克劳林公式给出：

${\displaystyle \ln(n!)-{\frac {\ln n}{2}}=\ln 1+\ln 2+\cdots +\ln(n-1)+{\frac {\ln n}{2}}=n\ln n-n+1+\sum _{k=2}^{m}{\frac {B_{k}{(-1)}^{k}}{k(k-1)}}\left({\frac {1}{n^{k-1}}}-1\right)+R_{m,n}.}$

其中B_k是伯努利数，R_m,n是欧拉-麦克劳林公式中的余项。取极限，可得：

${\displaystyle \lim _{n\to \infty }\left(\ln n!-n\ln n+n-{\frac {\ln n}{2}}\right)=1-\sum _{k=2}^{m}{\frac {B_{k}{(-1)}^{k}}{k(k-1)}}+\lim _{n\to \infty }R_{m,n}.}$

把这个极限记为y。由于欧拉-麦克劳林公式中的余项R_m,n满足：

${\displaystyle R_{m,n}=\lim _{n\to \infty }R_{m,n}+O\left({\frac {1}{n^{2m-1}}}\right),}$

其中用到了大O符号，与以上的方程结合，便得出对数形式的近似公式：

${\displaystyle \ln n!=n\ln \left({\frac {n}{e}}\right)+{\frac {\ln n}{2}}+y+\sum _{k=2}^{m}{\frac {B_{k}{(-1)}^{k}}{k(k-1)n^{k-1}}}+O\left({\frac {1}{n^{2m-1}}}\right).}$

两边取指数，并选择任何正整数m，便得到了一个含有未知数e^y的公式。当m=1时，公式为：

${\displaystyle n!=e^{y}{\sqrt {n}}~{\left({\frac {n}{e}}\right)}^{n}\left[1+O\left({\frac {1}{n}}\right)\right]}$

将上述表达式代入沃利斯乘积公式，并令n趋于无穷，便可以得出e^y（${\displaystyle e^{y}={\sqrt {2\pi }}}$）。因此，我们便得出斯特林公式：

${\displaystyle n!={\sqrt {2\pi n}}~{\left({\frac {n}{e}}\right)}^{n}\left[1+O\left({\frac {1}{n}}\right)\right]}$

这个公式也可以反复使用分部积分法来得出，首项可以通过最速下降法得到。把以下的和

${\displaystyle \ln(n!)=\sum _{j=1}^{n}\ln j}$

用积分近似代替，可以得出不含${\displaystyle {\sqrt {2\pi n}}}$的因子的斯特林公式（这个因子通常在实际应用中无关）：

${\displaystyle \sum _{j=1}^{n}\ln j\approx \int _{1}^{n}\ln x\,dx=n\ln n-n+1.}$

方法二

考虑如下近似关系：

${\displaystyle {\begin{aligned}\ln n!&=\sum _{k=1}^{n}\ln k\\&\approx \int _{1}^{n}\ln x\mathrm {d} x=n\ln n-n+1\\&\approx n\ln n-n\end{aligned}}}$

以及阶乘的积分形式（Γ函数）：

${\displaystyle n!=\int _{0}^{\infty }\mathrm {e} ^{-x}x^{n}\mathrm {d} x}$

其中被积函数取对数后求导结果可写成：

${\displaystyle {\frac {\mathrm {d} }{\mathrm {d} x}}\ln(\mathrm {e} ^{-x}x^{n})={\frac {\mathrm {d} }{\mathrm {d} x}}(n\ln x-x)={\frac {n}{x}}-1}$

因此被积函数只有在${\displaystyle x=n}$附近才趋于平稳，现令${\displaystyle x\equiv n+\xi }$且${\displaystyle \xi \ll n}$，并改写成如下形式：

${\displaystyle \ln(x^{n}\mathrm {e} ^{-x})=n\ln x-x=n\ln(n+\xi )-(n+\xi )}$

于是有：

${\displaystyle {\begin{aligned}\ln(n+\xi )&=\ln \left[n\left(1+{\frac {\xi }{n}}\right)\right]=\ln n+\ln \left(1+{\frac {\xi }{n}}\right)\\&=\ln n+{\frac {\xi }{n}}-{\frac {\xi ^{2}}{2n^{2}}}+\cdots \end{aligned}}}$

因此可以得到：

${\displaystyle {\begin{aligned}\ln(x^{n}\mathrm {e} ^{-x})&=n\ln(n+\xi )-(n+\xi )\\&=n\ln n+\xi -{\frac {\xi ^{2}}{2n}}-n-\xi +\cdots \\&=n\ln n-n-{\frac {\xi ^{2}}{2n}}+\cdots \end{aligned}}}$

两边同时作指数函数运算，变成：

${\displaystyle {\begin{aligned}x^{n}\mathrm {e} ^{-x}&\approx \mathrm {e} ^{n\ln n}\mathrm {e} ^{-n}\mathrm {e} ^{-{\frac {\xi ^{2}}{2n}}}\\&=n^{n}\mathrm {e} ^{-n}\mathrm {e} ^{-{\frac {\xi ^{2}}{2n}}}\end{aligned}}}$

代入到阶乘的积分表示形式，得到：

${\displaystyle {\begin{aligned}n!&\approx \int _{-n}^{\infty }n^{n}\mathrm {e} ^{-n}\mathrm {e} ^{-{\frac {\xi ^{2}}{2n}}}\mathrm {d} \xi \\&\approx n^{n}\mathrm {e} ^{-n}\int _{-\infty }^{\infty }\mathrm {e} ^{-{\frac {\xi ^{2}}{2n}}}\mathrm {d} \xi \end{aligned}}}$

约等号后侧得到的是高斯积分的形式，计算该积分得出斯特林公式：

${\displaystyle n!\approx n^{n}\mathrm {e} ^{-n}{\sqrt {2\pi n}}}$

收敛速率和误差估计

y轴表示截断的斯特林级数的相对误差，x轴表示所使用的项数。

更加精确的近似公式为：

${\displaystyle n!={\sqrt {2\pi n}}\;\left({\frac {n}{e}}\right)^{n}e^{\lambda _{n}}}$

其中：

${\displaystyle {\frac {1}{12n+1}}<\lambda _{n}<{\frac {1}{12n}}.}$

斯特林公式实际上是以下级数（现在称为斯特林级数）的第一个近似值：

${\displaystyle n!={\sqrt {2\pi n}}\left({n \over e}\right)^{n}\left(1+{1 \over 12n}+{1 \over 288n^{2}}-{139 \over 51840n^{3}}-{571 \over 2488320n^{4}}+\cdots \right).}$

当${\displaystyle n\to \infty }$时，截断级数的误差等于第一个省略掉的项。这是渐近展开式的一个例子。它不是一个收敛级数；对于任何特殊值n，级数的准确性只在取有限个项时达到最大，如果再取更多的项，则准确性将变得越来越差。

阶乘的对数的渐近展开式也称为斯特林级数：

${\displaystyle \ln n!=n\ln n-n+{1 \over 2}\ln(2\pi n)+{1 \over 12n}-{1 \over 360n^{3}}+{1 \over 1260n^{5}}-{1 \over 1680n^{7}}+\cdots .}$

在这种情况下，级数的误差总是与第一个省略掉的项异号，且最多同大小。

伽玛函数的斯特林公式

对于所有正整数，有：

${\displaystyle n!=\Pi (n)=\Gamma (n+1).}$

然而，伽玛函数与阶乘不一样，它对于所有复数都有定义。尽管如此，斯特林公式仍然适用。如果${\displaystyle \Re (z)>0}$，那么：

${\displaystyle \ln \Gamma (z)=\left(z-{\frac {1}{2}}\right)\ln z-z+{\frac {\ln {2\pi }}{2}}+2\int _{0}^{\infty }{\frac {\arctan {\frac {t}{z}}}{\exp(2\pi t)-1}}\,dt.}$

反复使用分部积分法，可得以下渐近展开式：

${\displaystyle \ln \Gamma (z)=\left(z-{\frac {1}{2}}\right)\ln z-z+{\frac {\ln {2\pi }}{2}}+\sum _{n=1}^{\infty }{\frac {(-1)^{n-1}B_{n}}{2n(2n-1)z^{2n-1}}}}$

其中B_n是第n个伯努利数。当${\displaystyle |\arg z|<\pi -\epsilon }$，其中ε是正数时，这个公式对于绝对值足够大的z是适用的，当使用了最初m个项时，误差项为${\displaystyle O\left(z^{-m-{\frac {1}{2}}}\right)}$。对应的近似值可以写为：

${\displaystyle \Gamma (z)={\sqrt {\frac {2\pi }{z}}}~{\left({\frac {z}{e}}\right)}^{z}\left[1+O\left({\frac {1}{z}}\right)\right].}$

斯特林公式的收敛形式

欲得出斯特林公式的一个收敛形式，我们必须计算：

${\displaystyle \int _{0}^{\infty }{\frac {2\arctan {\frac {t}{z}}}{\exp(2\pi t)-1}}\,dt=\ln \Gamma (z)-\left(z-{\frac {1}{2}}\right)\ln z+z-{\frac {1}{2}}\ln(2\pi ).}$

一种方法是利用含有上升阶乘幂的级数。如果${\displaystyle z^{\overline {n}}=z(z+1)\cdots (z+n-1)}$，那么：

${\displaystyle \int _{0}^{\infty }{\frac {2\arctan {\frac {t}{z}}}{\exp(2\pi t)-1}}\,dt=\sum _{n=1}^{\infty }{\frac {c_{n}}{(z+1)^{\overline {n}}}}}$

其中：

${\displaystyle c_{n}={\frac {1}{n}}\int _{0}^{1}x^{\overline {n}}\left(x-{\frac {1}{2}}\right)\,dx.}$

从中可以得出斯特林级数的一个收敛形式：

${\displaystyle \ln \Gamma (z)=\left(z-{\frac {1}{2}}\right)\ln z-z+{\frac {\ln {2\pi }}{2}}}$

${\displaystyle {}+{\frac {1}{12(z+1)}}+{\frac {1}{12(z+1)(z+2)}}+{\frac {59}{360(z+1)(z+2)(z+3)}}+{\frac {29}{60(z+1)(z+2)(z+3)(z+4)}}+\cdots }$

它在${\displaystyle \Re (z)>0}$时收敛。

适用于计算器的形式

以下的近似值

${\displaystyle \Gamma (z)\approx {\sqrt {\frac {2\pi }{z}}}\left({\frac {z}{e}}{\sqrt {z\sinh {\frac {1}{z}}+{\frac {1}{810z^{6}}}}}\right)^{z},}$

或

${\displaystyle 2\ln \Gamma (z)\approx \ln(2\pi )-\ln z+z\left[2\ln z+\ln \left(z\sinh {\frac {1}{z}}+{\frac {1}{810z^{6}}}\right)-2\right],}$

可以通过把斯特林公式整理，并注意到它的幂级数与双曲正弦函数的泰勒级数展开式的相似性来得出。当z的实数部分大于8时，这个近似值精确到小数点后8位。2002年，Robert H. Windschitl建议计算器用这个公式来计算伽玛函数。

Gergő Nemes在2007年提出了一个近似公式，它的精确度与Windschitl的公式相等，但更加简单：

${\displaystyle \Gamma (z)\approx {\sqrt {\frac {2\pi }{z}}}\left[{\frac {1}{e}}\left(z+{\frac {1}{12z-{\frac {1}{10z}}}}\right)\right]^{z},}$

或

${\displaystyle \ln \Gamma (z)\approx {\frac {1}{2}}\left[\ln(2\pi )-\ln z\right]+z\left[\ln \left(z+{\frac {1}{12z-{\frac {1}{10z}}}}\right)-1\right]}$

参考文献

• Abramowitz, M. and Stegun, I., Handbook of Mathematical Functions, http://www.math.hkbu.edu.hk/support/aands/toc.htm（页面存档备份，存于互联网档案馆）
• 夏慧异，章家顺, 史特灵公式在高等数学中的应用, http://ordostsg.org.cn:1080/KCMS/detail/detail.aspx?filename=XUSJ200904014&dbcode=CJFR&dbname=
• Paris, R. B., and Kaminsky, D., Asymptotics and the Mellin-Barnes Integrals, Cambridge University Press, 2001
• Whittaker, E. T., and Watson, G. N., A Course in Modern Analysis, fourth edition, Cambridge University Press, 1963. ISBN 0-521-58807-3
• Viktor T. Toth. "Programmable Calculators: Calculators and the Gamma Function". 2006 [2018-11-18]. （原始内容存档于2007-02-23）., modified 2006
• 埃里克·韦斯坦因. Stirling's Approximation. MathWorld.
• Stirling's approximation at PlanetMath.

1. Dutka, Jacques, The early history of the factorial function, Archive for History of Exact Sciences, 1991, 43 (3): 225–249, doi:10.1007/BF00389433
2. Le Cam, L., The central limit theorem around 1935, Statistical Science, 1986, 1 (1): 78–96, JSTOR 2245503, MR 0833276, doi:10.1214/ss/1177013818; see p. 81, "The result, obtained using a formula originally proved by de Moivre but now called Stirling's formula, occurs in his 'Doctrine of Chances' of 1733."
3. Pearson, Karl, Historical note on the origin of the normal curve of errors, Biometrika, 1924, 16 (3/4): 402–404 [p. 403], JSTOR 2331714, doi:10.2307/2331714, I consider that the fact that Stirling showed that De Moivre's arithmetical constant was ${\displaystyle {\sqrt {2\pi }}}$ does not entitle him to claim the theorem, [...]