Skip to content

Instantly share code, notes, and snippets.

@vvanirudh
Created April 14, 2024 22:05
Show Gist options
  • Save vvanirudh/c92a00bf482f174a5f1d203e1ba1a407 to your computer and use it in GitHub Desktop.
Save vvanirudh/c92a00bf482f174a5f1d203e1ba1a407 to your computer and use it in GitHub Desktop.
import matplotlib.pyplot as plt
import numpy as np
from scipy.stats import norm
import math
plt.rcParams["figure.figsize"] = (5,2)
fig, ax = plt.subplots()
pdf = lambda x: 0.7 * norm.pdf(x, 0.5, 0.1) + 0.3 * norm.pdf(x, 1, 0.1)
dx = 1e-2
xaxis = np.arange(0.25, 1.25, dx)
plt.plot(xaxis, list(map(pdf, xaxis)), label="p(x)", color='blue')
N = 100
np.random.seed(0)
def sample():
if np.random.rand() < 0.7:
return np.random.normal(0.5, 0.1)
return np.random.normal(1, 0.1)
samples = [sample() for _ in range(N)]
plt.plot(samples, np.zeros(len(samples)), 'x', label="samples", color='orange')
delta = 0.1
x = 1.0
ax.stem([x - delta/2, x + delta/2], [3, 3], markerfmt=' ', linefmt='k:', basefmt=' ')
k = lambda u: 1 if abs(u) <= 0.5 else 0
sorted_samples = sorted(samples)
plt.plot(sorted_samples, [k((x_n - x) / delta) for x_n in sorted_samples], label='k((x - x_n)/delta)', color='red')
K = lambda x: sum([k((x - x_n) / delta) for x_n in samples])
p = lambda x: K(x) / (N * delta)
plt.plot(xaxis, list(map(p, xaxis)), label="p(x) estimate", color='red')
box_at_x_n = lambda x, x_n: k((x - x_n) / delta)
sum_of_N_boxes_at_x = lambda x: sum([box_at_x_n(x, x_n) for x_n in samples])
p = lambda x: sum_of_N_boxes_at_x(x) / (N * delta)
plt.plot(xaxis, list(map(p, xaxis)), label="p(x) estimate", color='red')
x = 1.0
gaussian_k = lambda u: (1/(2*math.pi)**0.5) * math.exp(-0.5 * u**2)
plt.plot(sorted_samples, [gaussian_k((x_n - x) / delta) for x_n in sorted_samples], label='k((x - x_n)/delta)', color='red')
delta = 0.01
K = lambda x: sum([gaussian_k((x - x_n) / delta) for x_n in samples])
p = lambda x: K(x) / (N * delta)
plt.plot(xaxis, list(map(p, xaxis)), label="p(x) estimate", color='red')
plt.legend()
plt.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment