<jupyter_start><jupyter_text>Stirlingの公式（対数近似）* $\log n! \sim n\log n - n$* $n!$はおおよそ$\left(\frac{n}{e}\right)^n$になる* 参考: [スターリングの公式（対数近似）の導出](https://starpentagon.net/analytics/stirling_log_formula/)<jupyter_code>import numpy as np
import pandas as pd

import matplotlib.pyplot as plt<jupyter_output><empty_output><jupyter_text>$\log n!$の上からの評価<jupyter_code>MIN_X = 0.5
MAX_X = 10

x = np.linspace(MIN_X, MAX_X, 100)
y = np.log(x)

p = plt.plot(x, y, label='$\log x$')
p = plt.hlines([0], MIN_X, MAX_X)

p = plt.xlim(MIN_X, MAX_X-0.5)
p = plt.xticks(range(1, MAX_X+1))
p = plt.ylim([-0.2, 2.3])

# 面積log kの矩形を描画
for k in range(2, MAX_X):
    p = plt.vlines(k, 0, np.log(k), linestyles='dashed')
    p = plt.hlines(np.log(k), k, k+1, linestyles='dashed')

p = plt.legend()
plt.show(p)<jupyter_output><empty_output><jupyter_text>$\log n!$の下からの評価<jupyter_code>MIN_X = 0.5
MAX_X = 10

x = np.linspace(MIN_X, MAX_X, 100)
y = np.log(x)

p = plt.plot(x, y, label='$\log x$')
p = plt.hlines([0], MIN_X, MAX_X)

p = plt.xlim(MIN_X, MAX_X-0.5)
p = plt.xticks(range(1, MAX_X+1))
p = plt.ylim([-0.2, 2.3])

# 面積log kの矩形を描画
for k in range(2, MAX_X):
    p = plt.vlines(k-1, 0, np.log(k), linestyles='dashed')
    p = plt.hlines(np.log(k), k-1, k, linestyles='dashed')

p = plt.vlines(MAX_X-1, 0, np.log(MAX_X), linestyles='dashed')

p = plt.legend()
plt.show(p)<jupyter_output><empty_output><jupyter_text>$n \log n - n$の近似精度<jupyter_code>def log_factorial(n):
    '''log n!を返す'''
    val = 0.0
    
    for i in range(1, n+1):
        val += np.log(i)
        
    return val
# test of log_factorial
eps = 10**-5

assert abs(log_factorial(1) - 0.0) < eps
assert abs(log_factorial(2) - np.log(2)) < eps
assert abs(log_factorial(5) - np.log(120)) < eps
def log_factorial_approx(n):
    '''log n!の近似: n log n - nを返す'''
    return n * np.log(n) - n
# test of log_factorial_approx
assert abs(log_factorial_approx(1) - (-1)) < eps
assert abs(log_factorial_approx(2) - (2 * np.log(2) - 2)) < eps
# log_factorial, log_factorial_approxをplot
n_list = range(1, 50+1)
y_fact = [log_factorial(n) for n in n_list]
y_approx = [log_factorial_approx(n) for n in n_list]

p = plt.plot(n_list, y_fact, label='$\log n!$')
p = plt.plot(n_list, y_approx, label='$n \log n - n$')

p = plt.legend()

plt.show(p)
# 近似精度を評価
n_list = [5, 10, 20, 50, 100, 1000]

approx_df = pd.DataFrame()

approx_df['n'] = n_list
approx_df['log n!'] = [log_factorial(n) for n in n_list]
approx_df['n log(n)-n'] = [log_factorial_approx(n) for n in n_list]
approx_df['error(%)'] = 100 * (approx_df['log n!'] - approx_df['n log(n)-n']) / approx_df['log n!']
pd.options.display.float_format = '{:.1f}'.format
approx_df<jupyter_output><empty_output><jupyter_text>$n!$と$\left(\frac{n}{e}\right)^n$の比較<jupyter_code>n_list = [5, 10, 20, 50, 100]

approx_df = pd.DataFrame()

approx_df['n'] = n_list
approx_df['n!'] = [np.exp(log_factorial(n)) for n in n_list]
approx_df['(n/e)^n'] = [np.exp(log_factorial_approx(n)) for n in n_list]
approx_df['error(%)'] = 100 * (approx_df['n!'] - approx_df['(n/e)^n']) / approx_df['n!']
pd.options.display.float_format = None
pd.options.display.precision = 2

approx_df<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1977294228949858, 'fraction_numerical': 0.03122043519394513, 'mean_word_length': 4.705035971223022, 'pattern_counts': {'":': 0, '<': 26, '<?xml version=': 0, '>': 21, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '6244', 'license': 'MIT', 'n_tokens_mistral': 1602, 'n_tokens_neox': 1559, 'n_words': 281, 'path': 'math/Stirling/01_Stirling_log_formula.ipynb', 'repo_name': 'starpentagon/python_scripts', 'size': '3171'}