スターリングの近似

出典: フリー百科事典『ウィキペディア(Wikipedia)』
移動: 案内検索
(ln n!) と (n ln nn) は n が増大するにつれて近づいていく

スターリングの近似: Stirling's approximation)またはスターリングの公式: Stirling's formula)は、階乗、或いはガンマ関数漸近近似である。数学者ジェイムズ・スターリング英語版に因む。

応用上よく使われる形の公式は、ランダウの記号を用いて、

\log(n!) = n\log(n) - n +O(\log(n))

である。O(log(n)) における次の項は (1/2)log(2πn) である。よって公式のより正確な形は

n! \sim \sqrt{2\pi n} \left(\frac{n}{e}\right)^n

である。これは漸近公式英語版であり、スターリングの近似は

\lim_{n\to\infty} \frac{n!}{\sqrt{2\pi n} (n/e)^n} = 1

という性質を持つ。

n! の漸近近似よりもむしろ上下からの評価が必要なことがある。任意の正の整数 n に対して、

\sqrt{2\pi}\ n^{n+1/2}e^{-n} \le n! \le e\ n^{n+1/2}e^{-n}

が成り立ち、従って任意の n ≥ 1 に対して比 \frac{n!}{n^{n+1/2}e^{-n}} は常に \sqrt{2\pi} = 2.5066... A019727e = 2.71828... A001113 の間にある。

ガンマ関数 \Gamma(z) は正の整数 n に対して \Gamma(n)=(n-1)! を満たす。スターリングの近似はガンマ関数に対しても考えることができ

\Gamma(z+1)\approx\sqrt{2\pi z}\left(\frac{z}{e}\right)^z\qquad(|\arg z|<\pi,|z|\gg0)

が成り立つ。これは

\lim_{z\to\infty}\frac{\Gamma(z+1)}{\sqrt{2\pi z}(z/e)^z}=1\qquad(|\arg z|<\pi)

を意味する。但し、|\arg z|=\pi に近づくにつれ近似の誤差が大きくなる(極限の収束が遅くなる)ため、応用上は反射公式などを用いて |\arg z|\le\pi/2 程度に制限することが多い。

導出[編集]

初等的な方法による導出[編集]

スターリングの公式の導出を厳密にするためにはオイラーの和公式、或いは鞍点法など複素解析の技法を用いる必要があるが、多少の無理を許せば初等的に導くことも可能である。先ず、階乗の対数を積分で近似する。

\log(n!)=\sum_{k=1}^{n}{\log k}
 \approx\sum_{k=1}^{n}\int_{k-1/2}^{k+1/2}\log x\, dx
 =\int_{1/2}^{n+1/2}\log x\, dx

\int_{1/2}^{n+1/2}\log x\, dx=\left[x\log{x}-x\right]_{1/2}^{n+1/2}
 =\left(n+\frac{1}{2}\right) \log\left(n+\frac{1}{2}\right)-\left(n+\frac{1}{2}\right)
 -\left(\frac{1}{2}\right) \log\left(\frac{1}{2}\right)+\frac{1}{2}

然し、

 \begin{align}
 \log1-\int_{1/2}^{3/2}\log x\, dx &\simeq0.045 \\
 \log2-\int_{3/2}^{5/2}\log x\, dx &\simeq0.011 \\
 \log3-\int_{5/2}^{7/2}\log x\, dx &\simeq0.005 \\
\end{align}

などという誤差があるから、定数項と誤差をまとめて log c とし、c の値は後に求める。

\log(n!)\approx\log{c}+\left(n+\frac{1}{2}\right)\log\left(n+\frac{1}{2}\right)-\left(n+\frac{1}{2}\right)

両辺の真数を取ると

n!\approx c\left(n+\frac{1}{2}\right)^{n+1/2}\left(\frac{1}{e}\right)^{n+1/2}
 =ce^{-1/2}\left(n+\frac{1}{2}\right)^{1/2}\left(1+\frac{1}{2n}\right)^n \left(\frac{n}{e}\right)^n

となり

\lim_{n\to\infty}\left(n+\frac{1}{2}\right)^{1/2}=n^{1/2}=\sqrt{n}

\lim_{n\to\infty}\left(1+\frac{1}{2n}\right)^n=e^{1/2}

であるから

n!\approx c\sqrt{n}\left(\frac{n}{e}\right)^n

を得る。ウォリス積を変形して

\frac{\pi}{2}=\prod_{k=1}^{\infty}\frac{(2k)^2}{(2k+1)(2k-1)}\approx\frac{1}{2n+1}\left(\frac{(2n)!!}{(2n-1)!!}\right)^2=\frac{1}{2n+1}\left(\frac{(2n)!!(2n)!!}{(2n)!}\right)^2=\frac{2^{4n}(n!)^4}{(2n+1)((2n)!)^2}

これを代入すると

\frac{\pi}{2}=\frac{2^{4n}c^4n^2(n/e)^{4n}}{(2n+1)c^2(2n)(2n/e)^{4n}}=\frac{c^2n^2}{(2n+1)(2n)}\approx\frac{c^2}{4}

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

を得る。

精度を改善するために

n!\approx\sqrt{2\pi n}\left(\frac{n}{e}\right)^n\left(1+\frac{\alpha}{n}+\frac{\beta}{n^2}+\frac{\gamma}{n^3}\right)

とすると

\begin{align}(n+1)!
 &\approx\sqrt{2\pi{(n+1)}}\left(\frac{n+1}{e}\right)^{n+1}\left(
 1+\frac{\alpha}{n+1}+\frac{\beta}{(n+1)^2}+\frac{\gamma}{(n+1)^3}\right)\\
 &\approx\sqrt{2\pi{(n+1)}}\left(\frac{n+1}{e}\right)^{n+1}\left(
 1+\frac{\alpha}{n}-\frac{\alpha}{n^2}+\frac{\alpha}{n^3}-\frac{\alpha}{n^4}
 +\frac{\beta}{n^2}-\frac{2\beta}{n^3}+\frac{3\beta}{n^4}
 +\frac{\gamma}{n^3}-\frac{3\gamma}{n^3}+O\left(\frac{1}{n^5}\right)\right)
\end{align}

となる。然し

(n+1)!=(n+1)n!\approx(n+1)\sqrt{2\pi n}\left(\frac{n}{e}\right)^n
 \left(1+\frac{\alpha}{n}+\frac{\beta}{n^2}+\frac{\gamma}{n^3}\right)

であるから

\left(\frac{n+1}{n}\right)^{n+1/2}\frac{1}{e}
 \left(1+\frac{\alpha}{n}-\frac{\alpha}{n^2}+\frac{\alpha}{n^3}-\frac{\alpha}{n^4}
 +\frac{\beta}{n^2}-\frac{2\beta}{n^3}+\frac{3\beta}{n^4}
 +\frac{\gamma}{n^3}-\frac{3\gamma}{n^3}+O\left(\frac{1}{n^5}\right)\right)
 =1+\frac{\alpha}{n}+\frac{\beta}{n^2}+\frac{\gamma}{n^3}

となるはずである。再び対数を用いる。

\begin{align}\left(n+\frac{1}{2}\right)\log\left(\frac{n+1}{n}\right)-1
 &\approx\left(n+\frac{1}{2}\right)
 \left(\frac{1}{n}-\frac{1}{2n^2}+\frac{1}{3n^3}-\frac{1}{4n^4}+\frac{1}{5n^5}-O\left(\frac{1}{n^6}\right)\right)-1\\
 &\approx\frac{1}{12n^2}-\frac{1}{12n^3}+\frac{3}{40n^4}+O\left(\frac{1}{n^5}\right)
\end{align}

\left(\frac{n+1}{n}\right)^{n+1/2}\frac{1}{e}
\approx1+\frac{1}{12n^2}-\frac{1}{12n^3}+\frac{3}{40n^4}+\frac{1}{2(12n^2)^2}+O\left(\frac{1}{n^5}\right)

これを左辺に代入すると

1+\frac{\alpha}{n}+\frac{1}{n^2}\left( -\alpha+\beta+\frac{1}{12} \right)
 +\frac{1}{n^3}\left( \frac{13}{12}\alpha-2\beta+\gamma-\frac{1}{12} \right)
 +\frac{1}{n^4} \left( -\frac{7}{6}\alpha+\frac{37}{12}\beta-3\gamma+\frac{113}{1440} \right)
 +O\left(\frac{1}{n^5}\right)=1+\frac{\alpha}{n}+\frac{\beta}{n^2}+\frac{\gamma}{n^3}

となり、両辺の係数を比較して α=1/12,β=1/288,γ=-130/51840 を得る。結局、

n!\approx\sqrt{2\pi n}\left(\frac{n}{e}\right)^n\left(1+\frac{1}{12n}+\frac{1}{288n^2}-\frac{139}{51840n^3}\right)

を得る。

オイラーの和公式による導出[編集]

オイラーの乗積表示によるガンマ関数の定義の対数をとり

\begin{align}\log\Gamma(z)&=\log\big((z-1)\Gamma(z-1)\big)\\
&=\lim_{N\to\infty}(z-1)\log{N}+\sum_{n=1}^{N}\big(\log{n}-\log{(n+z-1)}\big)\\
\end{align}

f(n)=\log{n}-\log{(n+z-1)}オイラーの和公式を適用すれば

\begin{align}\log\Gamma(z)
&=\lim_{N\to\infty}(z-1)\log{N}+\int_{n=1}^{N}f(n)dn+\frac{1}{2}\big(f(N)+f(1)\big)\\
&\qquad+\sum_{k=1}^{m}\frac{B_{2k}}{(2k)!}\left(f^{(2k-1)}(N)-f^{(2k-1)}(1)\right)+\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)}{(2m+1)!}f^{(2m+1)}(n)dn\\
&=\lim_{N\to\infty}(z-1)\log{N}+\bigg[{n}\log{n}-{(n+z-1)}\log(n+z-1)\bigg]_{n=1}^{N}+\frac{1}{2}\big(\log{N}-\log{(N+z-1)}-\log{z}\big)\\
&\qquad+\sum_{k=1}^{m}\frac{B_{2k}}{(2k)(2k-1)}\left(\frac{1}{N^{2k-1}}-\frac{1}{(N+z-1)^{2k-1}}-1+\frac{1}{z^{2k-1}}\right)+\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)}{2m+1}\left(\frac{1}{n^{2m+1}}-\frac{1}{(n+z-1)^{2m+1}}\right)dn\\
&=\lim_{N\to\infty}\left(N+z-\frac{1}{2}\right)\big(\log{N}-\log{(N+z-1)}\big)+\left(z+\frac{1}{2}\right)\log{z}\\
&\qquad+\sum_{k=1}^{m}\frac{B_{2k}}{(2k)(2k-1)}\left(\frac{1}{N^{2k-1}}-\frac{1}{(N+z-1)^{2k-1}}-1+\frac{1}{z^{2k-1}}\right)+\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)}{2m+1}\left(\frac{1}{n^{2m+1}}-\frac{1}{(n+z-1)^{2m+1}}\right)dn\\
&=-z+1+\left(z-\frac{1}{2}\right)\log{z}-\sum_{k=1}^{m}\frac{B_{2k}}{(2k)(2k-1)}\left(1-\frac{1}{z^{2k-1}}\right)+\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)}{2m+1}\left(\frac{1}{n^{2m+1}}-\frac{1}{(n+z-1)^{2m+1}}\right)dn\\
\end{align}

となる。右辺の定数を集めて

C=1-\sum_{k=1}^{m}\frac{B_{2k}}{(2k)(2k-1)}+\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)dn}{(2m+1)n^{2m+1}}

とすれば

\log\Gamma(z)=C-z+\left(z-\frac{1}{2}\right)\log{z}+\sum_{k=1}^{m}\frac{B_{2k}}{(2k)(2k-1)z^{2k-1}}-\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)dn}{(2m+1)(n+z-1)^{2m+1}}

となり、この主要部をガンマ関数の反射公式に代入してz\to{i\infty}とすれば

\Gamma(z)\Gamma(1-z)=-z\Gamma(z)\Gamma(-z)=\frac{\pi}{\sin\pi{z}}
{\pi}i+\log{z}+\log\Gamma(z)+\log\Gamma(-z)-\log\sin{{\pi}z}-\log{\pi}=0
{\pi}i+\log{z}+C-z+\left(z-\frac{1}{2}\right)\log{z}+C+z+\left(-z-\frac{1}{2}\right)\left({\pi}i+\log{z}\right)-\log\sin{{\pi}z}-\log{\pi}=0
2C+\frac{{\pi}i}{2}-{\pi}iz-\log\sin{{\pi}z}-\log{\pi}=0

となるが

\log\sin{{\pi}z}=\log\left(e^{{\pi}iz}-e^{-{\pi}iz}\right)-\frac{{\pi}i}{2}-\log{2}\approx{\pi}iz-\frac{{\pi}i}{2}-\log{2}

であるから

C=\frac{\log{2\pi}}{2}

を得る。剰余項については

\alpha=\frac{2}{1+\cos(\arg{z})}

として

\begin{align}
\left|\int_{n=1}^{N}\frac{B_{2m+1}(n-\lfloor{n}\rfloor)dn}{(2m+1)(n+z-1)^{2m+1}}\right|
&\le\frac{\left|B_{2m}\right|}{2m}\int_{n=1}^{N}\frac{dn}{\left|n+z-1\right|^{2m+1}}\\
&\le\frac{\left|B_{2m}\right|}{2m\alpha^{2m+1}}\int_{n=1}^{N}\frac{dn}{(n+\left|z\right|-1)^{2m+1}}
=\frac{\left|B_{2m}\right|}{\alpha^{2m+1}\left|z\right|^{2m}}=O\left(z^{-2m}\right)\\
\end{align}

である。故に

\log\Gamma(z)=\frac{\log{2\pi}}{2}-z+\left(z-\frac{1}{2}\right)\log{z}+\sum_{k=1}^{m}\frac{B_{2k}}{(2k)(2k-1)z^{2k-1}}+O\left(z^{-2m}\right)

を得る。最初の数項を書き下せば

\log\Gamma(z)\approx\log\sqrt{2\pi}-z+\left(z-\frac{1}{2}\right)\log{z}
 +\frac{1}{12z}-\frac{1}{360z^3}+\frac{1}{1260z^5}-\frac{1}{1680z^7}+\frac{1}{1188z^9}

\Gamma(z)\approx\sqrt{\frac{2\pi}{z}}\left(\frac{z}{e}\right)^z
 \exp\left(\frac{1}{12z}-\frac{1}{360z^3}+\frac{1}{1260z^5}-\frac{1}{1680z^7}+\frac{1}{1188z^9}\right)

とやり、指数関数のテイラー展開により

\Gamma(z)\approx\sqrt{\frac{2\pi}{z}}\left(\frac{z}{e}\right)^z
 \left(1+\frac{1}{12z}+\frac{1}{288z^2}-\frac{139}{51840z^3}\right)

となる。

鞍点法による導出[編集]

スターリングの公式は鞍点法の好適例とされることが多いが、実際に複素平面全体(負の実数を除く)で漸近近似が成立することを鞍点法によって示すのは困難であるから、ここではzを正の実数に限定する。ガンマ関数t=z(1+u)の置換により

\begin{align}\Gamma(z+1)
&=\int_{0}^{\infty}{t^ze^{-t}dt}\\
&=\int_{-1}^{\infty}{z^z(1+u)^ze^{-z-zu}zdu}\\
&=z^{z+1}e^{-z}\int_{-1}^{\infty}{e^{-z\left(u-\log(1+u)\right)}du}\\
&=z^{z+1}e^{-z}\left(\int_{-1}^{-\epsilon}{e^{-z\left(u-\log(1+u)\right)}du}+\int_{-\epsilon}^{\epsilon}{e^{-z\left(u-\log(1+u)\right)}du}+\int_{\epsilon}^{\infty}{e^{-z\left(u-\log(1+u)\right)}du}\right)\qquad(\epsilon\ll1)\\
\end{align}

となるが、zが十分に大きければu=0の附近が支配的であるから

\Gamma(z+1)\approx{z^{z+1}e^{-z}\int_{-\epsilon}^{\epsilon}{e^{-z\left(u-\log(1+u)\right)}du}}\approx{z^{z+1}e^{-z}\int_{-\epsilon}^{\epsilon}{e^{-zu^2/2}du}}\approx{z^{z+1}e^{-z}\int_{-\infty}^{\infty}{e^{-zu^2/2}du}}

という近似が許され、ガウス積分により

\Gamma(z+1)\approx{z^{z+1}e^{-z}\sqrt\frac{2\pi}{z}=\sqrt{2{\pi}z}\left(\frac{z}{e}\right)^z}

を得る。\epsilon=z^{-1/3}として、近似の誤差は

\begin{align}\left|\int_{-1}^{-\epsilon}e^{-z\left(u-\log(1+u)\right)}du\right|
&\le\frac{1}{\epsilon}\int_{-1}^{-\epsilon}\left|\frac{-u}{(1+u)}e^{-z\left(u-\log(1+u)\right)}\right|du=\frac{1}{\epsilon{z}}\left[e^{-z\left(u-\log(1+u)\right)}\right]_{-1}^{-\epsilon}\\
&=\frac{1}{\epsilon{z}}\left(e^{-z\left(\epsilon^2/2+O(\epsilon^3)\right)}\right)-0\approx{z^{2/3}}e^{-z^{1/3}/2}
\end{align}
\begin{align}\left|\int_{\epsilon}^{\infty}e^{-z\left(u-\log(1+u)\right)}du\right|
&\le\frac{2}{\epsilon}\int_{\epsilon}^{\infty}\left|\frac{2u}{(1+u)\epsilon}e^{-z\left(u-\log(1+u)\right)}\right|du=\frac{2}{\epsilon{z}}\left[-e^{-z)\left(u-\log(1+u)\right)}\right]_{\epsilon}^{\infty}\\
&=0-\frac{2}{\epsilon{z}}\left(e^{-z)\left(\epsilon^2/2+O(\epsilon^3)\right)}\right)\approx{2z^{2/3}}e^{-z^{1/3}/2}
\end{align}
\begin{align}\int_{-\epsilon}^{\epsilon}e^{-z\left(u-\log(1+u)\right)}du
&=\int_{-\epsilon}^{\epsilon}e^{-z\left(u^2/2-u^3/3+O(u^4)\right)}du=\int_{-\epsilon}^{\epsilon}e^{-zu^2/2}e^{z\left(u^3/3+O(u^4)\right)}du\\
&=\int_{-\epsilon}^{\epsilon}e^{-zu^2/2}\left(1+zu^3/3+zO(u^4)\right)du\\
&=\int_{-\epsilon}^{\epsilon}e^{-zu^2/2}\left(1+zO(u^4)\right)du\\
\end{align}
\left|\int_{-\epsilon}^{\epsilon}e^{-z\left(u-\log(1+u)\right)}du-\int_{-\epsilon}^{\epsilon}e^{-zu^2/2}du\right|\le\left|zO(\epsilon^4)\int_{-\epsilon}^{\epsilon}e^{-zu^2/2}du\right|=\left|O(z^{-1/3})\int_{-\epsilon}^{\epsilon}e^{-zu^2/2}du\right|

であり

{z^{2/3}}e^{-z^{1/3}/2}\ll{z^{-1/2}}\qquad(z\to\infty)

であるから

\Gamma(z+1)=\sqrt{2{\pi}z}\left(\frac{z}{e}\right)^z\left(1+O(z^{-1/3})\right)

を得る。これは

\lim_{z\to\infty}\frac{\Gamma(z+1)}{\sqrt{2{\pi}z}\left(\frac{z}{e}\right)^z}=1\qquad(|\arg{z}|<{\pi})

を示すに十分である。但し、実際の誤差はO(z^{-1})であるが、それを鞍点法で示すのは困難である。

収束の速度と誤差見積もり[編集]

より正確に記すと、次のようになる。

n! = \sqrt{2 \pi n} \; \left(\frac{n}{e}\right)^{n}e^{\lambda_n}

ここで

\frac{1}{12n+1} < \lambda_n < \frac{1}{12n}.

スターリングの公式は以下の級数(スターリング級数)の近似(初項で打ち切ったもの)である。


  n!=\sqrt{2\pi n}\left({n\over e}\right)^n
  \left(
   1
   +{1\over12n}
   +{1\over288n^2}
   -{139\over51840n^3}
   -{571\over2488320n^4}
   + \cdots
  \right).

n \to \infty としたとき、省かれた級数はその最初の項とそれ以降が相殺するように漸近していく。これは漸近展開の一例である。

以下のような階乗の対数の漸近展開も「スターリング級数」と呼ぶ。


  \ln n!=n\ln n - n + {1\over 2}\ln(2\pi n)
   +{1\over12n}
   -{1\over360n^3}
   +{1\over1260n^5}
   -{1\over 1680n^7}
   +\cdots.

この場合、誤差は打ち切った級数の初項と同じ符号で同程度の大きさであることが知られている。

ガンマ関数に対するスターリングの公式[編集]

すべての正の整数に対して、

n! = \Pi(n) = \Gamma(n+1)

が成り立つ。ここで Γ はガンマ関数を表す。

しかしながら、パイ関数英語版 \Pi(z)=\Gamma(z+1) は、階乗とは異なり、より広く、正でない整数を除いてすべての複素数に対して定義される。それにもかかわらず、スターリングの公式をなお適用することができる。Re(z) > 0 であれば

\log(\Gamma (z)) =\left(z-\tfrac{1}{2}\right)\log(z) -z + \tfrac{1}{2}\log(2 \pi)
 +2\int_0^\infty \frac{\arctan(t/z)}{\exp(2\pi t)-1}\,{\rm d}t

部分積分を繰り返すことで次が得られる

\log(\Gamma (z)) \sim \left(z-\tfrac{1}{2}\right)\log(z) -z + \tfrac{1}{2}\log(2 \pi)
 +\sum_{n=1}^\infty \frac{B_{2n}}{2n(2n-1)z^{2n-1}}

ここで Bnn 番目のベルヌーイ数である。(無限和は収束しないので、この公式は漸近展開にすぎないことに注意する。)公式はεを正数として |arg(z)| < π−ε であるときに絶対値の十分大きい z に対して成り立つ。最初の m 項が使われるとき誤差項は O(z^{-2m-1}) である。対応する近似は

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

のように書ける。この漸近展開のより進んだ応用は Re(z) が定数の複素変数 z に対してである。例えば直線 1/4 + it 上でリーマン・ジーゲルテータ関数英語版の Im(z) において適用されたスターリングの公式を見よ。

ビネーの公式[編集]

スターリングの公式は収束しない級数を伴うので解析的に扱いづらいが、収束しない級数を収束する積分に換えたものとしてビネーの(第二)公式がある。

\Gamma(z+1)=\sqrt{2\pi z}\left(\frac{z}{e}\right)^ze^{\mu(z)},\quad
\mu(z)=2\int_0^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt,\quad(\real{z}>0)

ビネーの公式は、スターリングの級数を形式的に(収束条件を無視して)操作することによっても導かれるが、厳密には対数ガンマ関数の導関数にアベル・プラナの和公式を適用して得られる。

\log\Gamma(z)=\lim_{N\to\infty}z\log{N}+\sum_{n=1}^{N}\log{n}-\sum_{n=0}^{N}\log(n+z)

\frac{d}{dz}\log\Gamma(z)=\lim_{N\to\infty}\log{N}-\sum_{n=0}^{N}\frac{1}{n+z}

\frac{d^2}{dz^2}\log\Gamma(z)=\sum_{n=0}^{\infty}\frac{1}{(n+z)^2}

\real{z}>0ならf=(n+z)^{-2}は右半平面において正則であるからプラナの和公式により

\begin{align}\frac{d^2}{dz^2}\log\Gamma(z)
 &=\int_0^\infty \frac{1}{(n+z)^2}dt+\frac{1}{2z^2}+i\int_0^\infty \frac{(z+it)^{-2}-(z-it)^{-2}}{e^{2\pi t}-1}dt\\
 &=\frac{1}{z}+\frac{1}{2z^2}+i\int_0^\infty \frac{(z+it)^{-2}-(z-it)^{-2}}{e^{2\pi t}-1}dt\\
\end{align}

積分して

\frac{d}{dz}\log\Gamma(z)=C_1+\log{z}-\frac{1}{2z}
 +i\int_0^\infty \frac{-(z+it)^{-1}+(z-it)^{-1}}{e^{2\pi t}-1}dt

\begin{align}\log\Gamma(z)
 &=C_2+C_1z+z\log z-\frac{\log z}{2}+i\int_0^\infty \frac{-\log{(z+it)}+\log{(z-it)}}{e^{2\pi t}-1}dt\\
 &=C_2+C_1z+\left(z-\frac{1}{2}\right)\log z+2\int_0^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt\\
\end{align}

\real{z}>0なら|\tan(t/z)|<Mは有界であるから

\begin{align}
\left|\int_0^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt\right|
 &\le\left|\int_0^{|z|^{1/2}}\frac{\arctan(t/z)}{e^{2\pi t}-1}dt\right|
 +\left|\int_{|z|^{1/2}}^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt\right|\\
 &\le\int_0^{|z|^{1/2}}\frac{\sum_{k=1}^\infty |t/z|^kdt}{2\pi{t}}
 +M\int_{|z|^{1/2}}^\infty \frac{e^{2\pi}dt}{(e^{2\pi}-1)e^{t}e^{2\pi}}\qquad(|z|>1)\\
 &\le\int_0^{|z|^{1/2}}\frac{dt}{2\pi(|z|-t)}+\frac{Me^{2\pi}}{e^{2\pi}-1}\int_{t=|z|^{1/2}}^\infty e^{-2\pi t}dt\\
\end{align}

\lim_{z\to\infty}\left|\int_0^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt\right|
 \le\lim_{z\to\infty}\frac{1}{2\pi}\Big[\log\left(|z|-t\right)\Big]_0^{|z|^{1/2}}
 +\lim_{z\to\infty}\frac{Me^{2\pi}}{2\pi\left(e^{2\pi}-1\right)}\Big[-e^{-2\pi t}\Big]_{|z|^{1/2}}^\infty
=0

である。 スターリングの公式と比較して積分定数を求め

\log\Gamma(z)=\frac{1}{2}\log2\pi-z+\left(z-\frac{1}{2}\right)\log{z}
 +2\int_0^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt

真数に直して

\Gamma(z+1)=z\Gamma(z)=\sqrt{2\pi z}\left(\frac{z}{e}\right)^ze^{\mu(z)},\quad
\mu(z)=2\int_0^\infty \frac{\arctan(t/z)}{e^{2\pi t}-1}dt,\quad(\real{z}>0)

を得る。なお、ビネーの公式を元にして部分積分を繰り返すとスターリングの級数が得られる。

収束級数形式のスターリングの公式[編集]

トーマス・ベイズの John Canton への書簡が1763年王立協会により公表されている。それによると、スターリングの公式は収束級数ではないとされていた[1]

スターリングの公式の収束級数形式を得るには以下を評価する。

\int_0^\infty \frac{2\arctan(t/z)}{\exp(2\pi t)-1}\, dt
 =\ln\Gamma (z) - \left( z-\frac12 \right) \ln z +z - \frac12\ln(2\pi)

1つの方法として、階乗冪の逆数の収束級数を使う方法がある。z^{\overline n} = z(z+1) \cdots (z+n-1) としたとき、次のようになる。

\int_0^\infty \frac{2\arctan(t/z)}{\exp(2\pi t)-1} \, dt
 =\sum_{n=1}^\infty \frac{c_n}{(z+1)^{\overline n}}

ここで

c_n =\frac{1}{n} \int_0^1 x^{\overline n} \left( x-\frac{1}{2} \right)\, dx

である。以上から次のようなスターリング級数が得られる。

 \begin{align}
\ln \Gamma (z) &= \left( z-\frac{1}{2} \right) \ln z -z + \frac{1}{2}\ln(2\pi) \\
 &\quad +\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
\end{align}

これは、\Re(z)>0 のとき収束する。

計算機向けの変形[編集]

ガンマ関数の(関数電卓などの)計算機向けの近似として次の式がある。

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

これは、次と同等である。

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 がリソースの制限された計算機(電卓など)でのそれなりの正確性を持った近似としてこれを示した(参考文献の 'Toth' 参照)。

Gergő Nemes は 2007年にほぼ同程度の結果を与える近似式を提案した。こちらはより単純である。

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

これは、次と同等である。

 \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).

歴史[編集]

この公式を最初に発見したのはアブラーム・ド・モアブルで、次のような形式だった。

n!\sim [{\rm constant}]\cdot n^{n+1/2} e^{-n}

スターリングの貢献は定数が \sqrt{2\pi} であることを示したことである。より正確な形式は Jacques Binet が見出した。

スターリングの近似の「一次」バージョン n!=n^n は、マックス・プランクが1901年の黒体放射の論文で使用した。これは多量の光子や振動子についての黒体放射エネルギーの方程式にリンクしている。この近似は量子論でよく使われ、例えばピーター・デバイルイ・ド・ブロイも使っている。アルベルト・アインシュタインサティエンドラ・ボースは違う方式を採用した。非常に大きな n について確率分布をグラフに描画してみると、両者はほぼ平行になる。

外部リンク[編集]

参考文献[編集]