File:Donsker theorem for normal distributions.gif

import numpy as np import matplotlib.pyplot as plt import scipy

import tempfile import os import imageio

1. generate random samples from

def plot_bridge(samples=None, n=1000, dist=scipy.stats.norm, dist_name="normal", low=-3, high=+3):

   fig, axes = plt.subplot_mosaic("ABB", figsize=(20, 6))
   ax1 = axes["A"]
   ax2 = axes["B"]

   if samples is None:
       samples = dist.rvs(size=n)
   else:
       n = len(samples)
   x = sorted(samples)
   y = dist.cdf(x)
   ref_x = np.linspace(low,high, 1000)
   ref_y = dist.cdf(ref_x)

   # calculate empirical cdf
   ecdf = np.arange(n) / len(samples)
   # plot cdf of standard normal distribution and empirical cdf of samples
   ax1.plot(ref_x, ref_y, label='Standard CDF')
   ax1.plot(x, ecdf, label='Empirical CDF')
   ax1.set_title()
   ax1.set_xlim(low,high)
   ax1.legend()
   ax1.set_ylim(0,1)
   ax2.plot([0.0] + dist.cdf(x)                        .tolist() + [1.0], 
            [0.0] + (np.sqrt(n) * (ecdf - dist.cdf(x))).tolist() + [0.0])
   ax2.set_title('centered, scaled, and re-timed')
   ax2.set_ylim(-0.9, 0.9)
   fig.suptitle(f"{dist_name}, with n = {n}")
   return fig

def interpolate_counts(counts, frames_per_step):

   interpolated_counts = [counts[0]]
   for i in range(1,len(counts)):
       interval = (counts[i] - counts[i-1]) // i
       interpolated_counts += list(range(counts[i-1], counts[i], interval))
   return interpolated_counts + [counts[-1]]

with tempfile.TemporaryDirectory() as temp_dir:

   dist = scipy.stats.norm
   dist_name = "normal"
   low, high = -3.5, +3.5
   

   n_steps = 16
   frames_per_step = 10
   sample_counts = interpolate_counts([2**n for n in range(n_steps)], frames_per_step)
   n_frames = len(sample_counts)-1
   samples = dist.rvs(size=sample_counts[0]).tolist()
   
   for i in range(n_frames):
       samples += dist.rvs(size=sample_counts[i+1]-sample_counts[i]).tolist()
       fig = plot_bridge(samples, dist=dist, dist_name=dist_name, low=low,high=high)
       filename = os.path.join(temp_dir, f"plot_{i:03d}.png")
       fig.savefig(filename)
       plt.close(fig)

   # Compile images into GIF
   fps = 12
   images = []
   for i in range(n_frames):
       filename = os.path.join(temp_dir, f"plot_{i:03d}.png")
       images.append(imageio.imread(filename))
   imageio.mimsave(f"{dist_name} Donsker theorem.gif", images, duration=1/fps)