スターリングの近似

出典: フリー百科事典『ウィキペディア(Wikipedia)』
移動先: 案内検索
log n! と n log nnn → ∞ のとき漸近する

スターリングの近似: 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

である[1]。(ここで記号 \sim は両辺の比が(n → ∞ のとき) 1 に収束することを意味する。)

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

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

が成り立ち、従って任意の n に対して比 n!/(nn+1/2en)2π = 2.5066...e = 2.71828... の間にある。

スターリングの近似は階乗の複素引数への拡張の一つであるガンマ関数 Γ(z)(正の整数 n に対し Γ(n) = (n − 1)! が成り立つ;ボーア・モレルップの定理も参照)に拡張することができ、

\Gamma(z)\sim \sqrt{\frac{2\pi}{z}} \left( \frac{z}{e} \right)^z \qquad (\vert \arg z\vert \le \pi-\varepsilon ,\;\vert z\vert \to \infty )

が成り立つ(ただし ε > 0[2]|arg z | → π のときは収束が遅くなるため、応用上は相補公式などを用いて |arg z | ≤ π/2 程度に制限することが多い。

導出[編集]

初等的な導出[編集]

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

\log n!=\sum_{k=1}^n \log k\approx \sum_{k=1}^n \int_{k-\tfrac{1}{2}}^{k+\tfrac{1}{2}} \log x\, dx=\int_{\tfrac{1}{2}}^{n+\tfrac{1}{2}} \log x\, dx
\begin{align}
\int_{\tfrac{1}{2}}^{n+\tfrac{1}{2}} \log x\, dx
&=\Bigl[ x\log x-x\Bigr]_{\tfrac{1}{2}}^{n+\tfrac{1}{2}} \\
&=\left( n+\frac{1}{2} \right) \log \left( n+\frac{1}{2} \right) -\left( n+\frac{1}{2} \right)
 -\frac{1}{2} \log \frac{1}{2} +\frac{1}{2}
\end{align}

誤差

\begin{align}
\log 1-\int_{\tfrac{1}{2}}^{\tfrac{3}{2}} \log x\, dx &\approx 0.045 \\
\log2-\int_{\tfrac{3}{2}}^{\tfrac{5}{2}} \log x\, dx &\approx 0.011 \\
\log3-\int_{\tfrac{5}{2}}^{\tfrac{7}{2}} \log x\, dx &\approx 0.005 \\
\vdots &
\end{align}

があるから、誤差と定数項の和を log c とし、c の値は後に求める。

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

両辺の真数を取ると

\begin{align}
n! &=c\left( n+\frac{1}{2} \right)^{n+\tfrac{1}{2}} \left( \frac{1}{e} \right)^{n+\tfrac{1}{2}} \\
 &=ce^{-\tfrac{1}{2}} \left( n+\frac{1}{2} \right)^{\tfrac{1}{2}} \left( 1+\frac{1}{2n} \right)^n \left( \frac{n}{e} \right)^n \\
 &\sim c\sqrt{n} \left( \frac{n}{e} \right)^n
\end{align}

ウォリスの公式の系:

\frac{(2n)!}{n!^2} \sim \frac{4^n}{\sqrt{\pi n}}

に代入すると、

\frac{\sqrt{2} \cdot 4^n}{c\sqrt{n}} \sim \frac{4^n}{\sqrt{\pi n}}

\therefore c= \sqrt{2\pi}

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

を得る。

精度の改善[編集]

精度を改善するために

n!\sim \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)! &\sim \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\} \\
 &\sim\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!\sim (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+\tfrac{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
 &\sim \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 \\
 &\sim \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+\tfrac{1}{2}}\frac{1}{e}
\sim1+\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!\sim \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\} \\
 &\qquad+\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\} \\
 &\qquad +\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) \\
&\qquad+\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\sim \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)\sim \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)\sim \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)\sim \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^z e^{-t} dt \\
 &=\int_{-1}^{\infty} z^z (1+u)^z e^{-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}^{-\varepsilon} e^{-z\left\{ u-\log (1+u)\right\}}du+\int_{-\varepsilon}^{\varepsilon} {e^{-z\left\{ u-\log (1+u)\right\}} du}+\int_{\varepsilon}^{\infty} e^{-z \left\{ u-\log (1+u)\right\}} du \right) \qquad (\varepsilon \ll1)\\
\end{align}

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

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

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

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

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

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

であり

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

であるから

\Gamma (z+1)=\sqrt{2 \pi z} \left( \frac{z}{e} \right)^z \left( 1+O(z^{-\tfrac{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! \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)

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

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

\ln n!\sim n\ln n-n+\frac{1}{2} \ln 2\pi n+\frac{1}{12n} -\frac{1}{360n^3} +\frac{1}{1260n^5} -\frac{1}{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)}{e^{2\pi t}-1}\,dt

が成り立つ[3]。 部分積分を繰り返すことで次が得られる

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

ここで Bnn 番目のベルヌーイ数である。(無限和は収束しないので、この公式は漸近展開にすぎないことに注意する。)公式はεを正数として |arg z| < π − ε であるときに絶対値の十分大きい z に対して成り立つ。公式の右辺に現れる級数はスターリング級数と呼ばれる[4]。最初の 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 において適用されたスターリングの公式を見よ。

ビネーの公式[編集]

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

\Gamma (z+1)=\sqrt{2\pi z} \left( \frac{z}{e} \right)^z e^{\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_1 z+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_1 z+\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|^k dt}{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} \log 2\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)^z e^{\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年王立協会により公表されている。それによると、スターリングの公式は収束級数ではないとされていた[6]

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

\int_0^\infty \frac{2\arctan(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

一つの方法として、階乗冪の逆数の収束級数を使う方法がある。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)\sim \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)\sim \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 がリソースの制限された計算機(電卓など)でのそれなりの正確性を持った近似としてこれを示した[7]

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

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

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

\ln \Gamma (z) \sim \frac{1}{2} \left( \ln 2\pi -\ln z\right) +z\left\{ \ln \left( z+\frac{1}{12z-\frac{1}{10z}} \right) -1\right\}

歴史[編集]

この公式は最初に次の形でアブラーム・ド・モアブルにより1730年に発見された[8]

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

スターリングの貢献は定数が \sqrt{2\pi} であることを示したことである[8]。より正確な形式はジャック・ビネが見出した。

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

脚注[編集]

参考文献[編集]

外部リンク[編集]