The Gaussian Distribution

The Gaussian Distribution#

The Normal distribution#

The normal (or Gaussian) distribution is a ubiquitous one. It appears over and over again. There are two explanations as to why it appears so often:

It is the distribution of maximum uncertainty that matches a known expectation and a known variance variance.
It is the distribution that arises when you add a lot of random variables together.

We will learn about both these in subsequent lectures.

We write:

\[ X | \mu, \sigma \sim N(\mu, \sigma^2), \]

and we read “\(X\) conditioned on \(\mu\) and \(\sigma^2\) follows a normal distribution with expected value \(\mu\) and variance \(\sigma^2\).

Note

First, a lot of people write \(N(\mu, \sigma)\) instead of \(N(\mu, \sigma^2)\). You should always check what the author means.

Second, many people refer to \(\mu\) as the mean instead of the expected value. We will use the terms “mean” and “expected value” interchangeably when dealing with Gaussian distributions.

The standard normal distribution#

When we have zero mean and unit variance, we say that we have a standard normal distribution. If \(Z\) is a standard normal random variable, we write:

\[ Z\sim N(0,1). \]

The PDF of the standard normal is typically denoted by \(\phi(z)\) and is given by:

\[ \phi(z) = \frac{1}{\sqrt{2\pi}}\exp\left\{-\frac{z^2}{2}\right\}. \]

The CDF of the standard normal is denoted by \(\Phi(z)\) and is given by:

\[ \Phi(z) := p(Z \le z) = \int_{-\infty}^z \phi(z')dz'. \]

The CDF of the standard normal is not available in closed form. In the old days, we used to look it up in a table. Nowadays, we use a computer to calculate it.

Here is how you can get the PDF of the standard normal. First, let’s make a standard normal random variable in scipy.stats:

import scipy.stats as st
Z = st.norm()

You can evaluate it anywhere you want:

print(f"phi(0.5) = {Z.pdf(0.5):.2f}")

phi(0.5) = 0.35

Let’s plot the PDF:

../_images/6e12e45cc4070bfc0bf2d5eae67648cc8d52c012ff6459ba9abf76693e0dcdc1.svg

And here is the CDF:

../_images/84939d813874e36f2fc126d439349d9586257567bfab5caf9b89f934aa54e333.svg

Here is the expectation:

print(f"E[Z] = {Z.expect():.2f}")

E[Z] = -0.00

And the variance:

print(f"V[Z] = {Z.var():.2f}")

V[Z] = 1.00

Here is the probability that Z is between two numbers:

a = 1.0
b = 3.0
prob_Z_in_ab = Z.cdf(b) - Z.cdf(a)
print(f"p({a:.2f} <= Z <= {b:.2f}) = {prob_Z_in_ab:.2f}")

p(1.00 <= Z <= 3.00) = 0.16

And here is how you can sample:

Z.rvs(100)

And, of course, you can also sample using the functionality of numpy:

np.random.randn(100)

The general normal distribution#

Take

\[ X \sim N(\mu, \sigma^2). \]

The PDF is given by:

\[ f_X(x) = \frac{1}{\sqrt{2\pi\sigma^2}}\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}. \]

Note

Some people write \(N(x|\mu, \sigma^2)\) to refer to the PDF of \(N(\mu, \sigma^2)\) evaluated at \(x\).

If you mediate a little bit with this expression, you will notice that you can connect it to the standard normal PDF:

\[ f_X(x) = \frac{1}{\sigma}\phi\left(\frac{x-\mu}{\sigma}\right). \]

This highlights two things. First, the mean \(\mu\) shifts the distribution to the right or to the left.

# Plot the PDF of a standard normal and the PDF of a non-zero mean, unit variance normal
# with an arrow indicating the change from the first to the second.
fig, ax = plt.subplots()
zs = np.linspace(-4.0, 4.0, 100)
ax.plot(zs, Z.pdf(zs), label="$\phi(x)$")
ax.plot(zs, st.norm(loc=1.0).pdf(zs), label="$N(x|\mu=1,\sigma^2=1)$")
ax.set_xlabel("$x$")
ax.set_ylabel("$\phi(x)$")
ax.legend(loc="best", frameon=False)
ax.annotate("", xytext=(0.0, 0.4), xy=(1.0, 0.4), arrowprops=dict(arrowstyle="->",
            color=sns.color_palette()[2],
            lw=2))
plt.title("The mean translates the standard normal PDF")
sns.despine(trim=True);

../_images/7f114d0d8ec89ad37797906d095e41d4cc7115c2cbf6e8d4c98f23f8090c4fb9.svg

The standard deviation \(\sigma\) scales the distribution.

# Plot the PDF of a standard normal and the pdf of a zero mean Normal with variance two
fig, ax = plt.subplots()
zs = np.linspace(-4.0, 4.0, 100)
ax.plot(zs, Z.pdf(zs), label="$\phi(x)$")
ax.plot(zs, st.norm(loc=0.0, scale=np.sqrt(2.0)).pdf(zs), label="$N(x|\mu=0,\sigma^2=2)$")
ax.set_xlabel("$x$")
ax.set_ylabel("$p(x)$")
ax.legend(loc="best", frameon=False)
plt.title("The variance scales the standard normal PDF")
sns.despine(trim=True);

../_images/cf42dd148dd82f5d78bab816b9636e4b4d9232eb443c2d06ebce7b0c90856de4.svg

Quantiles of the normal#

There are a few more exciting things to know about the standard normal. For example, how can you find a value \(z_q\) such that the probability of \(Z\) being less than \(z_q\) is \(q\)? Mathematically, you wish to find this:

\[ \Phi(z_q) = q. \]

The point \(z_q\) is called the \(q\)-quantile. To find it, you need to do this:

\[ z_q = \Phi^{-1}(q). \]

For example, \(z_{0.50}\) is called the median (and it coincides with the expectation here). Another set of exciting quantiles is \(z_{0.025}\) and \(z_{0.975}\). Why? Because the probability that \(Z\) lies between them is \(95\)%. Here it is:

\[ p(z_{0.025} \le Z \le z_{0.975}) = \Phi(z_{0.975}) - \Phi(z_{0.025}) = 0.975 - 0.025 = 0.95. \]

Let’s find these quantiles and visualize them using the functionality of scipy.stats. We will use the percent point function (ppf), which the inverse of the CDF:

z_025 = Z.ppf(0.025) 
z_500 = Z.ppf(0.5)
z_975 = Z.ppf(0.975)
print(f"0.025 quantile of Z = {z_025:.2f}")
print(f"0.50 quantile of Z = {z_500:.2f}")
print(f"0.975 quantile of Z = {z_975:.2f}")

025 quantile of Z = -1.96
50 quantile of Z = 0.00
975 quantile of Z = 1.96

Here is how much probability there is between the two extreme quantiles:

print(f"p({z_025:.2f} <= Z <= {z_975:.2f}) = {Z.cdf(z_975) - Z.cdf(z_025):.2f}")

p(-1.96 <= Z <= 1.96) = 0.95

Note

Observe that the quantiles are symmetric around the mean (which is zero here). The 0.025 quantile is at -1.96 and the 0.975 quantile is at 1.96. Being engineers, we are going to call the -1.96 a two.

Let’s also visualize the quantiles on top of the PDF:

../_images/cab05fbac5392caf21845008542df68cfe91a6708ed38bc329eb1d3d1f672336.svg

Note

Quantiles vs percentiles A percentile is a quantile expressed as a percentage. So, you say the 0.025 quantile is the 2.5 percentile. The 0.975 quantile is the 97.5 percentile. And so on.

Questions#

Modify the code above so that you find and vizualize \(z_{0.001}\) and \(z_{0.999}\).
What is the difference between \(z_{0.999}\) and \(z_{0.001}\)?
What is the probability that \(Z\) is between \(z_{0.001}\) and \(z_{0.999}\)?

Getting any normal from the standard normal#

Knowledge of the quantiles of the standard normal is sufficient to give us the quantiles of any normal. Let \(Z\) be a standard normal. Take some \(\mu\) and \(\sigma\). Then, the random variable

\[ X = \mu + \sigma Z, \]

follows a \(N(\mu,\sigma^2)\).

To show this mathematically, we need to show that the PDF of \(X\) is the one we expect.

Proof

First, write down the CDF of \(X\):

\[ F_X(x) = p(X \le x) = p(\mu + \sigma Z \le x) = p\left(Z \le \frac{x-\mu}{\sigma}\right) = \Phi\left(\frac{x-\mu}{\sigma}\right). \]

Now take the derivative of the CDF to get the PDF:

\[ f_X(x) = F'_X(x) = \frac{1}{\sigma}\phi\left(\frac{x-\mu}{\sigma}\right), \]

which is exactly what we wanted.

The formula is extremely useful. For example, you can use it to make samples from any normal distribution using samples from the standard normal. Here is how:

mu = 1.0
sigma = 0.1
X = st.norm(mu, sigma)
xs = np.linspace(mu - 6.0 * sigma, mu + 6.0 * sigma, 100)
x_samples = mu + sigma * Z.rvs(size=10000)

Compare the histogram of the samples with the PDF of the normal:

fig, ax = plt.subplots()
ax.hist(x_samples, density=True, alpha=0.5, label="Samples from $X = \mu + \sigma Z$")
ax.plot(xs, X.pdf(xs), label="PDF of $N(\mu, \sigma^2)$")
ax.set_xlabel("$x$")
ax.set_ylabel("$p(x)$")
plt.legend(loc="best", frameon=False)
sns.despine(trim=True);

../_images/dcddf11e64100ca30541f61bfbe1b8c0089da7eae67a7098b1e2875aefe03e98.svg

How can you find the quantiles of this normal? Well, you can simply use the functionality of scipy.stats. As an example, let’s find \(x_{0.025}\):

x_025 = X.ppf(0.025)
print(f"0.025-quantile of N({mu:.2f}, {sigma:.2f}^2) = {x_025:1.2f}")

0.025-quantile of N(1.00, 0.10^2) = 0.80

But we can also find this quantile by exploiting the connection between \(X\) and \(Z\). The definition of a quantile of \(X\) is:

\[ p(X \le x_q) = q. \]

But, since \(X=\mu+\sigma Z\), this is equivalent to:

\[ p(\mu + \sigma Z \le x_q) = q, \]

which becomes:

\[ p(\sigma Z \le x_q-\mu) = q, \]

and then:

\[ p\left(Z \le \frac{x_q-\mu}{\sigma}\right) = q. \]

This is just:

\[ \Phi\left(\frac{x_q-\mu}{\sigma}\right) = q, \]

which tells us that \(\frac{x_q-\mu}{\sigma}\) is the \(q\)-quantile of \(Z\), i.e.,

\[ z_q = \frac{x_q-\mu}{\sigma}. \]

Solving for \(x_q\), we get:

\[ x_q = \mu + \sigma z_q. \]

Let’s do a sanity check:

z_025 = Z.ppf(0.025)
print(f"mu + sigma * z_025 = {mu + sigma * z_025:.2f}")

mu + sigma * z_025 = 0.80

which is the same as what we found before. So, let’s find also the 0.975-quantile:

z_975 = Z.ppf(0.975)
x_975 = mu + sigma * z_975
print(f"0.975-quantile of N({mu:.2f}, {sigma:.2f}^2) = {x_975:1.2f}")

0.975-quantile of N(1.00, 0.10^2) = 1.20

Let’s visualize the quantiles like we did before:

../_images/109829a5b1a83a4eca92ecafd8576a10107c2243620c4295dd346bb61bb0cf50.svg

Now, let’s find the distance between \(x_{2.5}\) and \(x_{97.5}\) in terms of the standard deviation \(\sigma\). We have:

\[ x_{97.5} - x_{2.5} = \mu + \sigma z_{97.5} - \mu - \sigma z_{2.5} = \sigma (z_{97.5} - z_{2.5}). \]

This is:

print(f"x_975 - x_025 ~= sigma * {z_975 - z_025:.2f}")

x_975 - x_025 ~= sigma * 3.92

Okay. So we see that 95% of the probability is contained within a \(3.92\sigma\) interval. This interval is centered at the median (which here is the same as the mode and the expectation of the probability density). The number 3.92 could be nicer, so we will round up to 4 intervals. A 4\(\sigma\) interval about the mean gives us a bit more than 95% of the probability, but it’s simpler to remember.

Note

Remember For a normal random variable \(N(\mu,\sigma^2)\), the 95% probability interval is about \(\mu \pm 4\sigma\). In other words:

\[ p(\mu - 2\sigma < X < \mu + 2 \sigma) \approx 0.95, \]

Questions#

Write code that finds exactly how much probability there is between \(\mu - 2\sigma\) and \(\mu + 2\sigma\), i.e., find \(p(\mu - 2\sigma < X < \mu + 2 \sigma)\).
Modify the code you just written, to find how much probability there is in \(\mu - 3\sigma\) and \(\mu + 3\sigma\), i.e., find \(p(\mu - 3\sigma < X < \mu + 3 \sigma)\). This is six-sigmas interval about the mean. Have you ever heard of the six-sigma process improvement technique?

The Gaussian Distribution

Contents

The Gaussian Distribution#

The Normal distribution#

The standard normal distribution#

The general normal distribution#

Quantiles of the normal#

Questions#

Getting any normal from the standard normal#

Questions#