Math — Multivariate Probability

A progressive implementation of multivariate statistics — from mean vectors and covariance matrices to correlation and the multivariate normal distribution, building toward the Gaussian foundations of machine learning with NumPy.

---

Learning Objectives

---

Data Conventions

This module uses two matrix layouts — a key learning point about data conventions in ML:

The transpose flips between conventions. Pay attention to which axis you compute the mean along and how the covariance matrix multiplication is ordered.

---

Task-by-Task Reference

Each task below highlights the unique challenge it posed and the new technique introduced — techniques from earlier tasks are not repeated.

---

Task 0 — Mean & Covariance (0-mean_cov.py)

Challenge: Compute both the mean vector and covariance matrix of a multivariate dataset $X$ — implementing the sample covariance formula from scratch without using np.cov.

Approach: Given $X$ of shape $(n, d)$, compute the mean with np.mean(X, axis=0, keepdims=True) to get shape $(1, d)$. Center the data: $X - \mu$. The sample covariance is $\frac{1}{n-1} (X - \mu)^\top (X - \mu)$ — the matrix multiplication of the centered data's transpose with itself, scaled by Bessel's correction. Inputs are validated in order: 2D ndarray check → $n \ge 2$ check.

New techniques introduced:

Key takeaway: The covariance matrix $\Sigma$ captures how each pair of dimensions varies together — $\Sigma_{ij}$ is the covariance between dimension $i$ and dimension $j$. Bessel's correction ($n-1$) gives an unbiased estimate of the population covariance. Without keepdims=True, the mean collapses to 1D, breaking the broadcasting in the centering step.

---

Task 1 — Correlation Matrix (1-correlation.py)

Challenge: Convert a covariance matrix into a correlation matrix — normalizing covariances to the $[-1, 1]$ range so relationships between variables measured in different units become directly comparable.

Approach: Extract the diagonal variances $\sigma_i^2$ with np.diagonal(). Compute scaling factors $\sigma_i^{-1/2}$ (element-wise power). Build a diagonal matrix $D$ from these factors with np.diagflat(). The correlation matrix is the sandwich product $D \cdot C \cdot D$, which scales each covariance $\sigma_{ij}$ by $\frac{1}{\sqrt{\sigma_{ii} \cdot \sigma_{jj}}}$.

New techniques introduced:

Key takeaway: Correlation is standardized covariance — it strips away the units so you can compare relationships across variables measured on different scales (e.g., height in cm vs. weight in kg). The diagonal scaling trick $DCD$ is a clean one-liner that transforms any covariance matrix into its correlation form.

---

Task (class) — Multivariate Normal Distribution (multinormal.py)

Challenge: Implement a class that fits a multivariate normal distribution from data and computes the PDF for any point — combining mean, covariance, determinant, and inverse into the full multivariate Gaussian formula. This is the culmination of all multivariate probability concepts.

Approach: The constructor accepts data of shape $(d, n)$ (features × samples — note the transposed convention from Task 0). It computes the mean along axis 1 (per-feature) with keepdims=True for shape $(d, 1)$. The covariance is $\frac{1}{n-1} (X - \mu)(X - \mu)^\top$ — note the reversed multiplication order vs. Task 0 because of the transposed layout. The pdf(x) method computes the full multivariate normal density:

$f(x) = (2\pi)^{-d/2} \, |\Sigma|^{-1/2} \, \exp\!\left(-\frac{1}{2}(x - \mu)^\top \Sigma^{-1} (x - \mu)\right)$

using np.linalg.det() for the determinant, np.linalg.inv() for the precision matrix, and .item() to extract the scalar result from the 1×1 matrix.

New techniques introduced:

Key takeaway: The multivariate normal PDF generalizes the 1D Gaussian bell curve to $d$ dimensions. The covariance matrix $\Sigma$ controls the shape (spread) and orientation (correlation) of the bell. The Mahalanobis distance $(x - \mu)^\top \Sigma^{-1} (x - \mu)$ inside the exponent measures distance accounting for correlations — points on the same density contour have the same Mahalanobis distance. The determinant $|\Sigma|$ in the normalization scales the peak height: larger determinant → more spread → lower peak.

---

Technique Inventory

---

Multivariate Probability Pipeline

Data X ──→ Mean μ ──→ Covariance Σ ──→ Correlation P ──→ MultiNormal(μ, Σ) ──→ PDF(x)
  (n×d)    (1×d)        (d×d)            (d×d)              class               scalar

1. Mean locates the center of the data in $d$-dimensional space
2. Covariance captures the spread and pairwise relationships
3. Correlation standardizes to $[-1, 1]$ for interpretability
4. MultiNormal packages mean and covariance into a full probability distribution
5. PDF answers: "how likely is this point under the fitted distribution?"

---

Resources

• NumPy mean — Official Docs
• NumPy linalg.det — Official Docs
• NumPy linalg.inv — Official Docs
• NumPy diagonal — Official Docs
• Covariance Matrix — Wikipedia
• Correlation — Wikipedia
• Multivariate Normal Distribution — Wikipedia
• Mahalanobis Distance — Wikipedia