Covariance Matrix

Compute the covariance matrix of a dataset without np.cov.

For X of shape (n, d) where rows are samples and columns are features,

C[i, j] = sum over samples of (x_i - mean_i) * (x_j - mean_j) / (n - 1)

which is X_centered.T @ X_centered / (n - 1).

The n - 1 is Bessel's correction: dividing by n underestimates the spread, because the deviations are measured from the sample mean rather than the true mean. This is the matrix PCA decomposes, so its diagonal is the per feature variance and the off-diagonal entries say which features move together.

Return a (d, d) array. If n < 2 the covariance is undefined, so return a (d, d) array of zeros.

Function signature

def covariance_matrix(X):

Example

Input

X = np.array([[1.0, 2.0], [3.0, 6.0], [5.0, 10.0]])
covariance_matrix(X)

Output

array([[ 4.,  8.],
       [ 8., 16.]])

Reasoning

The column means are 3 and 6. The deviations are (-2, 0, 2) and (-4, 0, 4). Dividing the sums of products by n - 1 = 2 gives variances 4 and 16 and a covariance of 8.

Constraints and notes

• Divide by n - 1, not n.
• When n < 2, return a (d, d) array of zeros.
• Do not use np.cov.