Multivariate Gaussians

Multivariate Gaussians are often presented as intimidating objects. The density function looks kind of scary:

$$f(\mathbf{x}) = \frac{1}{\sqrt{(2\pi)^n \det(\boldsymbol{\Sigma})}} \exp\left( -\frac{1}{2} (\mathbf{x} - \boldsymbol{\mu})^{\top} \boldsymbol{\Sigma}^{-1} (\mathbf{x} - \boldsymbol{\mu}) \right).$$

But there’s a clearer perspective: any multivariate Gaussian can be disentangled into independent 1-dimensional Gaussians. More precisely, there’s always a basis in which a multivariate Gaussian has independent coordinates. This basis is given by the eigenbasis of the covariance matrix, and the eigenvalues correspond to the variances along each axis.

Random vectors

Let $\mathbf{x} = (x_1, \ldots, x_k)^{\top} \in \mathbb{R}^{k}$ be a random vector. The expectation and covariance of $\mathbf{x}$ can be written as follows:

$$\mathbb{E}[\mathbf{x}] = \left(\mathbb{E}[x_1], \ldots, \mathbb{E}[x_k]\right)^{\top},$$ $$\text{Cov}(\mathbf{x}) = \mathbb{E}\left[\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)^{\top}\right].$$

Note that the covariance matrix $\text{Cov}(\mathbf{x})$ is symmetric and positive semi-definite. To see this, let’s fix some vector $\mathbf{a} \in \mathbb{R}^{k}$ and consider $\mathbf{a}^{\top} \text{Cov}(\mathbf{x}) \mathbf{a}$:

$$\begin{aligned} \mathbf{a}^{\top} \text{Cov}(\mathbf{x}) \mathbf{a} &= \mathbf{a}^{\top} \mathbb{E}\left[\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)^{\top}\right] \mathbf{a} \\ &= \mathbb{E}\left[\mathbf{a}^{\top} \left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)^{\top} \mathbf{a}\right] \\ &= \mathbb{E}\left[(\mathbf{a}^{\top} (\mathbf{x} - \mathbb{E}[\mathbf{x}]))^2\right] \\ &\geq 0. \end{aligned}$$

If we multiply a random vector $\mathbf{x} \in \mathbb{R}^{k}$ by a matrix $\mathbf{A} \in \mathbb{R}^{n \times k}$, then the covariance of the resulting random vector $\mathbf{y} = \mathbf{A}\mathbf{x} \in \mathbb{R}^{n}$ is given by:

$$\begin{aligned} \text{Cov}(\mathbf{y}) &= \mathbb{E}\left[\left(\mathbf{y} - \mathbb{E}[\mathbf{y}]\right)\left(\mathbf{y} - \mathbb{E}[\mathbf{y}]\right)^{\top}\right] \\ &= \mathbb{E}\left[\left(\mathbf{A}\mathbf{x} - \mathbf{A}\mathbb{E}[\mathbf{x}]\right)\left(\mathbf{A}\mathbf{x} - \mathbf{A}\mathbb{E}[\mathbf{x}]\right)^{\top}\right] \\ &= \mathbb{E}\left[\mathbf{A}\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)^{\top} \mathbf{A}^{\top}\right] \\ &= \mathbf{A} \mathbb{E}\left[\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)\left(\mathbf{x} - \mathbb{E}[\mathbf{x}]\right)^{\top}\right] \mathbf{A}^{\top} \\ &= \mathbf{A} \text{Cov}(\mathbf{x}) \mathbf{A}^{\top}. \end{aligned}$$

Multivariate Gaussian distribution

Now we can define the multivariate Gaussian distribution.

Let $\mathbf{g} = (g_1, \ldots, g_k)^{\top} \in \mathbb{R}^{k}$ be a vector of i.i.d. standard normal random variables. Given a matrix $\mathbf{A} \in \mathbb{R}^{n \times k}$, the covariance of $\mathbf{A}\mathbf{g}$ is $\mathbf{A} \text{Cov}(\mathbf{g}) \mathbf{A}^{\top} = \mathbf{A} \mathbf{I} \mathbf{A}^{\top} = \mathbf{A} \mathbf{A}^{\top}$.

The distribution of $\mathbf{A}\mathbf{g}$ is called a multivariate normal distribution with covariance $\boldsymbol{\Sigma} = \mathbf{A} \mathbf{A}^{\top}$, denoted $\mathcal{N}(\mathbf{0}, \boldsymbol{\Sigma})$.

Since $\boldsymbol{\Sigma}$ is symmetric positive semi-definite, we can write its eigendecomposition as

$$\boldsymbol{\Sigma} = \mathbf{Q} \boldsymbol{\Lambda} \mathbf{Q}^{\top},$$

where $\mathbf{Q} \in \mathbb{R}^{n \times n}$ is an orthogonal matrix, and $\boldsymbol{\Lambda} \in \mathbb{R}^{n \times n}$ is a diagonal matrix with eigenvalues $\lambda_1, \ldots, \lambda_n$ on the diagonal.

Now we can take

$$\mathbf{A} = \mathbf{Q} \boldsymbol{\Lambda}^{1/2} \quad \text{or} \quad \mathbf{A} = \mathbf{Q} \boldsymbol{\Lambda}^{1/2} \mathbf{Q}^{\top}$$

so that $\mathbf{A} \mathbf{A}^{\top} = \mathbf{Q} \boldsymbol{\Lambda}^{1/2} \boldsymbol{\Lambda}^{1/2} \mathbf{Q}^{\top} = \mathbf{Q} \boldsymbol{\Lambda} \mathbf{Q}^{\top} = \boldsymbol{\Sigma}$.

For example, take $\mathbf{x} = \mathbf{Q} \boldsymbol{\Lambda}^{1/2} \mathbf{g}$ for $\mathbf{g} \sim \mathcal{N}(\mathbf{0}, \mathbf{I})$. Then $\mathbf{x} \sim \mathcal{N}(\mathbf{0}, \boldsymbol{\Sigma})$.

Visualizing the transformation $\mathbf{x} = \mathbf{Q}\boldsymbol{\Lambda}^{1/2}\mathbf{g}$

1. Start: g ~ N(0, I)

i.i.d. standard normals

→

2. Scale: Λ^1/2g

stretch by λ₁^1/2, λ₂^1/2

→

3. Rotate: QΛ^1/2g

align to q₁, q₂

λ₁ = 2.5

λ₂ = 0.4

θ = 30°

Disentangling the multivariate Gaussian

Now label $\mathbf{q}_1, \ldots, \mathbf{q}_n$ as the columns of $\mathbf{Q}$.

Expanding $\mathbf{x} = \mathbf{Q} \boldsymbol{\Lambda}^{1/2} \mathbf{g}$:

$$\mathbf{x} = \sqrt{\lambda_1} g_1 \mathbf{q}_1 + \ldots + \sqrt{\lambda_n} g_n \mathbf{q}_n = \sum_{i=1}^{n} \sqrt{\lambda_i} g_i \mathbf{q}_i.$$

In the orthonormal basis $\mathbf{q}_1, \ldots, \mathbf{q}_n$, the vector $\mathbf{x}$ has coordinates $(\sqrt{\lambda_1} g_1, \ldots, \sqrt{\lambda_n} g_n)$. These coordinates are independent 1-dimensional Gaussian random variables, with variances $\lambda_1, \ldots, \lambda_n$.

If $\boldsymbol{\Sigma}$ is singular (i.e. $\det \boldsymbol{\Sigma} = 0$), some eigenvalues will be zero, say $\lambda_{k+1} = \ldots = \lambda_n = 0$. The random vector is then concentrated on the subspace spanned by $\mathbf{q}_1, \ldots, \mathbf{q}_k$, and has no density on $\mathbb{R}^n$. On this subspace, $\mathbf{x}$ has density (in the eigenbasis coordinates):

$$f(x_1, \ldots, x_k) = \prod_{i=1}^{k} \frac{1}{\sqrt{2\pi\lambda_i}} \exp\left( -\frac{x_i^2}{2\lambda_i} \right).$$

References

References are listed in alphabetical order.

1. Lecture 4: Multivariate normal distribution and multivariate CLT [link]
   Dmitry Panchenko. 2006.