ガンマ関数 │ Python 数値計算入門

1. ガンマ関数
2. ディガンマ関数
3. ポリガンマ関数
4. 不完全ガンマ関数

ガンマ関数

　整数 $n$ について階乗 $n!$ は
　
\[n!=\begin{cases}1 & (n=0)\\[6pt]n(n-1)(n-2)\ \cdots\ 2\cdot 1 & (n\geq 1)\end{cases}\]
によって定義されますが、$n$ を実部が正となる複素数 $z$ にまで拡大定義した連続関数を ガンマ関数 とよびます。
　
\[\Gamma(z)=\int_{0}^{\infty} e^{-t}t^{z-1}\tag{1}dt\]
　$z$ が整数 $n$ であるとき、ガンマ関数と階乗の間には次のような関係があります。
　
\[\Gamma(n+1)=n!\tag{2}\]
　$z$ が非整数であっても、$\Gamma(z)$ は以下の公式によって再帰計算することができます。
　
\[\Gamma(z+1)=z\Gamma(z)\tag{3}\]
　また、$z=1/2$ におけるガンマ関数の値を (1) によって計算すると
　
\[\Gamma \left( \frac{1}{2}\right)=\sqrt{\pi}\tag{4}\]
となり、半整数値におけるガンマ関数の値は
　
\[\Gamma \left( n+\frac{1}{2}\right) =\frac{(2n-1)!!}{2^n}\sqrt{\pi} \tag{5}\]
で表せることが知られています。
　
　math.gamma(x) を使って、正数 x におけるガンマ関数の値を計算できます。

# GAMMA_FUNCTION

# In[1]

from math import gamma

# Γ(4)=3!
val = gamma(4)

print(val)
# 6.0

　Matplotlib を使って、ガンマ関数をグラフに描いてみます。ただし、NumPy にはガンマ関数を計算する関数が用意されていないので、frompyfunc() を使って、math.gamma() をユニバーサル関数に変換しておきます。

# In[2]

import numpy as np
import matplotlib.pyplot as plt

# x座標データを作成
x = np.linspace(0.01, 5, 257)

# math.gamma()をユニバーサル化
gamma_u = np.frompyfunc(gamma, 1, 1)

# y座標データを作成
y = gamma_u(x)

# ガンマ関数のグラフを描画
fig = plt.figure(figsize = (6, 6))
ax = fig.add_subplot(111)
ax.set_title("Gamma Function", size = 15)
ax.grid()
ax.set_xlim(0, 5)
ax.set_ylim(0, 10)
ax.set_xlabel("x", size = 14, labelpad = 10)
ax.set_ylabel("Gamma(x)", size = 14, labelpad = 10)
ax.plot(x, y, color = "darkblue")

python math.gamma()
　
　math.lgamma(x) はガンマ関数の自然対数を返します。

# In[3]

# log{Γ(3.5)}
val = lgamma(3.5)

print("{:3.f}".format(val))
# 1.201

　SciPy の special.gamma(z) は任意の複素数 z を受け取ってガンマ関数の値を返します。

# In[4]

from scipy.special import gamma

# Γ(2+3i)
val = gamma(2 + 3j)

print("{:.3f}".format(val))
# -0.082+0.092j

　Matplotlib を使って、全実数領域におけるガンマ関数の値をグラフにプロットしてみます。

# In[5]

# x座標データを作成
x = np.linspace(-5, 5, 2251)

# y座標データを作成
y = gamma(x)

# ガンマ関数のグラフを描画
fig = plt.figure(figsize = (6, 6))
ax = fig.add_subplot(111)
ax.set_title("Gamma Function", size = 15)
ax.grid()
ax.set_xlim(-5, 5)
ax.set_ylim(-10, 10)
ax.set_xlabel("x", size = 14, labelpad = 10)
ax.set_ylabel("Gamma(x)", size = 14, labelpad = 10)
ax.plot(x, y, color = "darkblue")
plt.show()

scipy.gamma() plot
　
　scipy.special.gammaln() はガンマ関数の自然対数を高精度で計算します。

# In[6]

from scipy.special import gammaln

# logΓ(10)
val = gammaln(10)

print("{:.3f}".format(val))
# 12.802

ディガンマ関数

　ガンマ関数の対数微分を ディガンマ関数 とよびます。
　
\[\psi(z)=\frac{d}{dz}\ln\Gamma(z)\tag{6}\]
　SciPy にはディガンマ関数を計算する scipy.special.psi(x) が用意されています。

# DIGAMMA_FUNCTION

# In[1]

from scipy.special import psi

# ψ(-1.5),ψ(1),ψ(2)
val = psi(-1.5)

print("{:.3f}".format(val))
# 0.423

　$\gamma=-\psi(1)$ はオイラーの定数とよばれます。

# In[2]

# オイラーの定数
g = -psi(1)

print("γ = {}".format(g))
# γ = 0.5772156649015329

　オイラーの定数は numpy.euler_gamma で呼び出すこともできます。

# In[3]

import numpy as np

# NumPyのオイラー定数
g = np.euler_gamma

print("γ = {}".format(g))
# 0.5772156649015329

　ディガンマ関数には
　
\[\psi(x+1)=\psi(x)+\frac{1}{x}\tag{7}\]　
という性質があるので、整数点におけるディガンマ関数は有限調和級数からオイラーの定数を差し引いた形で表されます。
　
\[\psi(n+1)=1+\frac{1}{2}+\frac{1}{3}+\ \cdots \ +\frac{1}{n}-\gamma\tag{8}\]
　実数範囲でディガンマ関数をプロットするコードを載せておきます。

# In[4]

import matplotlib.pyplot as plt

# x座標データを作成
x = np.linspace(-5, 5, 2251)

# y座標データを作成
y = psi(x)

# ディガンマ関数のグラフを描画
fig = plt.figure(figsize = (6, 6))
ax = fig.add_subplot(111)
ax.set_title("Digamma Function", size = 15)
ax.grid()
ax.set_xlim(-5, 5)
ax.set_ylim(-10, 10)
ax.set_xlabel("x", size = 14, labelpad = 10)
ax.set_ylabel("ψ(x)", size = 14, labelpad = 10)
ax.plot(x, y, color = "darkblue")

python digamma function plot

ポリガンマ関数

　ディガンマ関数 $\psi(z)$ の $1$ 階微分 $\psi^{(1)}(z)$、$2$ 階微分 $\psi^{(2)}(z)$、$3$ 階微分 $\psi^{(3)}(z)$ ... をそれぞれトリガンマ関数、テトラガンマ関数、ペンタガンマ関数 ... とよび、これらの総称を ポリガンマ関数 $\psi^{(n)}(z)$ と定義します。
　
\[\psi^{(n)}=\frac{d^{n+1}}{dz^{n+1}}\ln\Gamma(z)\tag{9}\]
　SciPy にはポリガンマ関数を計算する scipy.special.polygamma(n,x) が用意されています。

# POLYGAMMA__FUNCTION

# In[1]

from scipy.special import polygamma

# z=1におけるディガンマ関数の値 
di = polygamma(0, 1)

# z=1におけるトリガンマ関数の値 
tri = polygamma(1, 1)

# z=1におけるテトラガンマ関数の値 
tetra = polygamma(2, 1)

print("ψ(1): {:.3f}".format(a))
print("ψ1(1): {:.3f}".format(b))
print("ψ2(1): {:.3f}".format(c))

# ψ(1): -0.577
# ψ1(1): 1.645
# ψ2(1): -2.404

不完全ガンマ関数

　ガンマ関数の積分区間 $[0,\infty]$ を 2 つに分割して、不完全ガンマ関数を次のように定義します。
　
\[\begin{align*}\gamma(a,x)&=\int_0^xe^{-t}t^{a-1}\tag{10}\\[6pt]\Gamma(a,x)&=\int_x^{\infty}e^{-t}t^{a-1}=\Gamma(a)-\gamma(a,x)\tag{11}\end{align*}\]
　(10) を第1種不完全ガンマ関数 (lower incomplete gamma function)、(11) を第2種不完全ガンマ関数 (upper incomplete gamma function) とよびます。
　scipy.special.gammainc(a,x) は第1種不完全ガンマ関数 $\gamma(a,x)$ を計算します。

INCOMPLETE_GAMMA_FUNCTION

# In[1]

import numpy as np
from scipy.special import gammainc
import matplotlib.pyplot as plt

# 区間[0,12]を128分割
x = np.linspace(0, 12, 129)

# FigureとAxes
fig = plt.figure(figsize=(6.5, 5))
ax = fig.add_subplot(111)
ax.grid()
ax.set_title("Lower incomplete gamma function", fontsize=16)
ax.set_xlabel("x", fontsize=16)
ax.set_ylabel("γ(a, x)", fontsize=16)
ax.set_xlim(0, 12)
#ax.set_ylim(0, 10)

# 第1種不完全ガンマ関数をプロット
for a in range(1, 6):
    ax.plot(x, gammainc(a, x), label="a={}".format(a))

# 凡例を表示
ax.legend(loc="lower right")

第1種不完全ガンマ関数
　
　scipy.special.gammaincc(a,x) は第2種不完全ガンマ関数 $\Gamma(a,x)$ を計算します。

# In[2:

# 区間[0,12]を128分割
x = np.linspace(0, 12, 129)

# FigureとAxes
fig = plt.figure(figsize=(6.5, 5))
ax = fig.add_subplot(111)
ax.grid()
ax.set_title("Upper incomplete gamma function", fontsize=16)
ax.set_xlabel("x", fontsize=16)
ax.set_ylabel("Γ(a, x)", fontsize=16)
ax.set_xlim(0, 12)
#ax.set_ylim(0, 10)

# 第2種不完全ガンマ関数をプロット
for a in range(1, 6):
    ax.plot(x, gammaincc(a, x), label="a={}".format(a))

# 凡例を表示
ax.legend(loc="upper right")

第2種不完全ガンマ関数

コメント（1件）

あとりえこばとより:

2023年5月19日 5:22 PM

人工知能にガンマ関数の用途を尋ねました

　
　ChatGPTにガンマ関数について、あれこれ質問してみました。
　
[Q] ガンマ関数はどのような用途で用いられますか？

[A] ガンマ関数は、数学や統計学などの分野で広く使用される特殊関数です。以下に、ガンマ関数の一部の主な用途を示します。
1. 階乗の一般化: ガンマ関数は、非整数における階乗の一般化と見なすことができます。正の整数nに対してn! = (n-1)! × n という関係が成り立つように、ガンマ関数ではガンマ(z+1) = z × ガンマ(z) という性質を持ちます。
2. 確率分布: ガンマ関数は、確率分布の形状パラメータや尺度パラメータの表現に使用されます。例えば、ガンマ分布やベータ分布などの確率分布は、ガンマ関数を使用して定義されます。
3. 統計解析: ガンマ関数は、最尤推定やベイズ統計などの統計解析手法で使用されます。ガンマ関数は、確率密度関数や累積分布関数の正規化に使用されることがあります。
4. 数値計算: ガンマ関数は、数値計算や科学技術計算において広く使用されます。特に、ガンマ関数の値を計算するためのアルゴリズムや数表が利用可能であり、複雑な数式や積分計算において役立ちます。
5. 物理学や工学: ガンマ関数は、物理学や工学の多くの分野で使用されます。例えば、波動関数や確率振幅などの表現にガンマ関数が登場することがあります。
　これらはガンマ関数の一部の用途ですが、その他にも様々な応用があります。ガンマ関数は数学的な表現力と実用性を持ち、多くの数学者、統計学者、科学者によって幅広く利用されています。

[Q] ガンマ関数は誰が発見しましたか？
　
[A] ガンマ関数は、18世紀の数学者レオンハルト・オイラーによって初めて研究されました。オイラーは、複数の数学的な関数や数列に関する業績を残し、ガンマ関数の性質や特殊値についても研究しました。オイラーは1730年代からガンマ関数に関する研究を始め、1737年に「De progressionibus transcendentibus seu quarum termini generales algebraice dari nequeunt」という論文で初めてガンマ関数について記述しました。彼はガンマ関数を階乗の一般化として扱い、ガンマ関数の性質や積分表現について論じました。その後、オイラーの研究を基に、他の数学者たちもガンマ関数について研究を進めました。特にカール・フリードリッヒ・ガウスやアドリアン＝マリ・ルジャンドルなどの数学者たちは、ガンマ関数の性質や応用について重要な貢献をしました。したがって、ガンマ関数の発見と初期の研究はレオンハルト・オイラーによって行われ、その後の数学者たちの研究によって洗練されてきました。

返信