Having difficult interpreting the eigenvectors for a simple 3x2 matrix

Question

I calculated the eigenvectors and eigenvalues from a covariance matrix given a data matrix of 3 columns and 2 rows.

I am trying to interpret results but I can't understand on how to interpret them.

Create a 2x3 matrix:

# Create a 2x3 matrix

data = np.around(np.random.uniform(size=(2,3)) * 100)

The data looks as follows:

    [
      [ 4., 65., 77.],
      [68., 12., 89.]
    ]

# Here each row represents one data point 
# and columns represent the features in the data set
# So there are 3 features and 2 data points

Calculate the mean for each feature in the data set.

mean = np.mean(data, axis = 0)

Center the data around origin, by subtracting mean from the data set.

difference = np.subtract(data, mean)

Now, calculate the covariance matrix:

cov = np.dot(difference.T, difference)

The cov matrix looks as follows:

[
    [ 2048. , -1696. ,   384. ],
    [-1696. ,  1404.5,  -318. ],
    [  384. ,  -318. ,    72. ]
]

As I understand about the covariance matrix, it explains the variance between all feature-pairs. Since there are 3 features, it gives out a 3x3 matrix explaining the variance between all possible pairs.

Finally, calculate the eigenvectors and eigenvalues:

val, vec = np.linalg.eigh(cov)

The vec matrix looks as follows:

[
    [ 0.60999981,  0.21639063,  0.76228297],
    [ 0.77451164,  0.040441  , -0.63126559],
    [ 0.16742745, -0.97546892,  0.14292806]
]

How do I interpret the the vector matrix? I understand what are eigenvectors physically. They do not change in position when an object undergoes a transformation but only a scalar change by their eigenvalues.

What are some possible ways, I could use this vec matrix?

Answer 1

The eigenvectors can be used to transform your data into a coordinate system in which no covariance is there. Assume we have a $p$ dimensional multivariate normal distribution with

$$\mathcal{N}(\boldsymbol{x}|\boldsymbol{\mu}=\boldsymbol{0},\boldsymbol{\Sigma})=\dfrac{1}{\sqrt{(2\pi)^p\det\boldsymbol{\Sigma}}}\exp\left[-\dfrac{1}{2}\boldsymbol{x}^T\boldsymbol{\Sigma}\boldsymbol{x} \right].\quad (*)$$

As the covariance matrix is real and symmetric we know it is diagonalizable and that we can scale the eigenvectors to represent an orthonormal basis by the set of all eigenvectors. The eigenvalue equation for the covariance matrix is given by

$$\boldsymbol{\Sigma}\boldsymbol{v}_i=\lambda_i\boldsymbol{v}_i \quad \forall i=1,...,p.$$

We can combine all equations into a single matrix equation

$$\boldsymbol{\Sigma}[\boldsymbol{v}_1,\ldots,\boldsymbol{v}_p]=[\boldsymbol{v}_1,\ldots,\boldsymbol{v}_p]\text{diag}\left[\lambda_1,\ldots,\lambda_p \right].$$

If we call $\boldsymbol{V}=[\boldsymbol{v}_1,\ldots,\boldsymbol{v}_p]$ and $\boldsymbol{\Lambda}=\left[\lambda_1,\ldots,\lambda_p \right].$ With these definitions in hand we can write the eigenvalue equation as

$$\boldsymbol{\Sigma}\boldsymbol{V}=\boldsymbol{V}\boldsymbol{\Lambda}$$ $$\implies \boldsymbol{\Sigma}=\boldsymbol{V}\boldsymbol{\Lambda}\boldsymbol{V}^{-1}.$$

as $\boldsymbol{V}$ is orthogonal (consists of orthonormal vectors) we can write $\boldsymbol{V}^{-1}=\boldsymbol{V}^T.$ This implies

$$\boldsymbol{\Sigma}=\boldsymbol{V}\boldsymbol{\Lambda}\boldsymbol{V}^{T}.$$

We plug this into $(*)$ and introduce the new variable $\boldsymbol{z}=\boldsymbol{V}^T\boldsymbol{x}$.

$$\mathcal{N}(\boldsymbol{x}|\boldsymbol{\mu}=\boldsymbol{0},\boldsymbol{\Sigma})=\dfrac{1}{\sqrt{(2\pi)^p\det\boldsymbol{V}\boldsymbol{\Lambda}\boldsymbol{V}^T}}\exp\left[-\dfrac{1}{2}\boldsymbol{x}^T\boldsymbol{V}\boldsymbol{\Lambda}\boldsymbol{V}^{T}\boldsymbol{x} \right]$$ $$=\dfrac{1}{\sqrt{(2\pi)^p\det\boldsymbol{V}\boldsymbol{\Lambda}\boldsymbol{V}^T}}\exp\left[-\dfrac{1}{2}\boldsymbol{z}\boldsymbol{\Lambda}\boldsymbol{z} \right]$$ $$\implies \mathcal{N}(\boldsymbol{z}|\boldsymbol{\mu}=\boldsymbol{0},\boldsymbol{\Lambda})=\dfrac{1}{\sqrt{(2\pi)^p\det\boldsymbol{\Lambda}}}\exp\left[-\dfrac{1}{2}\boldsymbol{z}\boldsymbol{\Lambda}\boldsymbol{z} \right].$$

In the last step I used $\det \boldsymbol{ABC}=\det\boldsymbol{A}\det\boldsymbol{B}\det\boldsymbol{C}$, $\det\boldsymbol{A}=\det\boldsymbol{A}^T$ for orthogonal matrices and $|\det\boldsymbol{A}|=1$ for orthogonal matrices. Hence, we proved that we can use the eigenvectors to linearly transform our variables to obtain a new coordinate system which has a diagonal covariance matrix $\boldsymbol{\Lambda}$ (which is diagonal). This implies that we do not have any covariance anymore.

Remark 1: This diagonalization is what the principal components analysis is also doing under the hood.

Remark 2: There is one suboptimal part of your code. You should explicitly determine the sample_size = data.shape[0] and then calculate cov = 1 / (sample_size - 1) * np.dot(difference.T, difference).