import numpy as np
import matplotlib.pyplot as plt

def uniform_sums(size):
    def sum_of_12():
        return sum(np.random.uniform(size=12)) - 6
    
    return [sum_of_12() for i in range(size)]
    
num_draws = 100000

uniforms = uniform_sums(size=num_draws)
normal = np.random.normal(size=num_draws)

bins = np.linspace(-6, 6)
plt.hist(uniforms, label='sums of 12', bins=bins, histtype='step')
plt.hist(normal, label='normal', bins=bins, histtype='step')
plt.legend()

import scipy.integrate
import scipy.stats

def inverse_transform(x):
    def integrand(t):
        return 1 / (2 * np.pi) * (2 * np.sin(t/2) / t)**12 * np.cos(t*x)
    
    p, _ = scipy.integrate.quad(integrand, -np.infty, np.infty)
    return p

sum_of_uniform_pdf = np.vectorize(inverse_transform)
norm_pdf = scipy.stats.norm().pdf

xs, step= np.linspace(-6, 6, 500, retstep=True)
ps_uniform = sum_of_uniform_pdf(xs)
ps_normal = norm_pdf(xs)

ps_uniform.sum() * step

plt.plot(xs, ps_uniform, label="sum of uniforms")
plt.plot(xs, ps_normal, label="N(0, 1)")
plt.legend()
plt.title("Both pdfs")

plt.plot(xs, ps_uniform - ps_normal)
plt.title("difference")

plt.plot(xs, (ps_uniform - ps_normal) / ps_normal)
plt.title("relative difference")

def how_many_to_exceed(threshold):
    total = 0
    count = 0
    while total < threshold:
        total += np.random.uniform()
        count += 1
    return count

num_draws = 100000
draws = [how_many_to_exceed(1) for i in range(num_draws)]
running_average = np.cumsum(draws) / (np.arange(num_draws) + 1.)

plt.plot(running_average, label='Average draws')
plt.axhline(np.e, color='red', linestyle='--', label='e = {:.8f}'.format(np.e))
plt.ylim(1, 4)
plt.legend()
plt.gca().set_xscale('log')