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

# 1
Normal = st.norm(loc = 50_000, scale = 2_000)
Normal.cdf(48_000) # cumulative distribution function <=
Normal.ppf(1 - .9)

47436.8968689108

# 2
Normal = st.norm(loc = 10, scale = 2)
Normal.cdf(12) - Normal.cdf(8)# a
1 - Normal.cdf(14)

0.02275013194817921

N = st.norm()
x = np.linspace(-5, 5, 101)
fx = N.pdf(x)
plt.plot(x, fx)

fig = plt.figure()
ax = fig.add_subplot(111, projection="3d")
X, Y = np.meshgrid(x, x)
XY = np.dstack((X, Y))
MvN = st.multivariate_normal(np.zeros(2), np.diag(np.ones(2)))
Z = MvN.pdf(XY)
ax.plot_surface(X, Y, Z)

MvN.pdf(XY[1, 4])

np.float64(2.473314011763375e-11)

XY[1, 4]

array([-4.6, -4.9])

N.pdf(-4.6) * N.pdf(-4.9)

np.float64(2.4733140117633728e-11)