def test_cdf(self):
"""Test that the implementation of the gaussian value through the
cumulative distribution function works as expected.
"""
from scipy.stats import norm
from scipy.special import erf
start = -5
stop = 5
n_points = 9
centers = np.array([0])
sigma = 1
area_cum = []
area_pdf = []
# Calculate errors for dfferent number of points
for n_points in range(2, 10):
axis = np.linspace(start, stop, n_points)
# Calculate with cumulative function
dx = (stop - start)/(n_points-1)
x = np.linspace(start-dx/2, stop+dx/2, n_points+1)
pos = x[np.newaxis, :] - centers[:, np.newaxis]
y = 1/2*(1 + erf(pos/(sigma*np.sqrt(2))))
f = np.sum(y, axis=0)
f_rolled = np.roll(f, -1)
pdf_cum = (f_rolled - f)[0:-1]/dx
# Calculate with probability function
dist2 = axis[np.newaxis, :] - centers[:, np.newaxis]
dist2 *= dist2
f = np.sum(np.exp(-dist2/(2*sigma**2)), axis=0)
f *= 1/math.sqrt(2*sigma**2*math.pi)
pdf_pdf = f
true_axis = np.linspace(start, stop, 200)
pdf_true = norm.pdf(true_axis, centers[0], sigma) # + norm.pdf(true_axis, centers[1], sigma)
# Calculate differences
sum_cum = np.sum(0.5*dx*(pdf_cum[:-1]+pdf_cum[1:]))
sum_pdf = np.sum(0.5*dx*(pdf_pdf[:-1]+pdf_pdf[1:]))
area_cum.append(sum_cum)
area_pdf.append(sum_pdf)
# mpl.plot(axis, pdf_pdf, linestyle=":", linewidth=3, color="r")
# mpl.plot(axis, pdf_cum, linewidth=1, color="g")
# mpl.plot(true_axis, pdf_true, linestyle="--", color="b")
# mpl.show()
mpl.plot(area_cum, linestyle=":", linewidth=3, color="r")
mpl.plot(area_pdf, linewidth=1, color="g")
# mpl.plot(true_axis, pdf_true, linestyle="--", color="b")
mpl.show()
评论列表
文章目录