s = np.random.normal(size=N)
p, x = np.histogram(s, bins=n)
x = x[:-1] + (x[1] - x[0]) / 2
f = UnivariateSpline(x, p, s=n)
plt.plot(x, f(x))
plt.show()
Plotting probability density function by sample with matplotlib