0
$\begingroup$

I recently learned about different methods of PCA. I decided to manually implement PCA in Python with Eigendecomposition of cov(X) and the Singular Value Decomposition of X and compare the results. I heard that the main difference is that SVD should give a more accurate result while taking longer to execute. Here is my implementation:

import numpy as np
from scipy.linalg import svd,eig
SIZE=(3,3) #the data shape

def svd_(data):
    mean=np.mean(data,axis=0)
    centered=data-mean #center data using column means
    print(f'{centered=}')
    U, s, Vh=svd(centered)
    print(f'{U=}\n{s=}\n,{Vh=}')

def eig_(data):
    mean=np.mean(data,axis=0)
    centered=data-mean #center data using column means
    print(f'{centered=}')
    cov=np.cov(centered.T)
    eigen_values, eigen_vectors=eig(cov)
    print(f'{eigen_values=}\n{eigen_vectors=}')


def main():
    np.random.seed(42) #fixed seed
    data=np.random.random(size=SIZE)
    print(f'data is\n{data}')
    svd_(data)
    print()
    eig_(data)

if __name__=='__main__':
    main()

I get the following output:

data is
[[0.37454012 0.95071431 0.73199394]
 [0.59865848 0.15601864 0.15599452]
 [0.05808361 0.86617615 0.60111501]]

#svd_ below
centered=array([[ 0.03077938,  0.29307794,  0.23562612],
       [ 0.25489775, -0.50161772, -0.3403733 ],
       [-0.28567713,  0.20853978,  0.10474719]])
U=array([[ 0.41479721,  0.7032851 ,  0.57735027],
       [-0.81646137,  0.00758237,  0.57735027],
       [ 0.40166416, -0.71086748,  0.57735027]])
s=array([8.05432409e-01, 2.49318619e-01, 2.61896384e-17])
,Vh=array([[-0.38500217,  0.76341895,  0.5186182 ],
       [ 0.90910975,  0.21687009,  0.35564986],
       [ 0.15903707,  0.60840683, -0.77752707]])  

#eig_ below
centered=array([[ 0.03077938,  0.29307794,  0.23562612],
       [ 0.25489775, -0.50161772, -0.3403733 ],
       [-0.28567713,  0.20853978,  0.10474719]])
eigen_values=array([ 3.24360683e-01+0.j,  3.10798868e-02+0.j, -4.80648600e-18+0.j])
eigen_vectors=array([[ 0.38500217, -0.90910975, -0.15903707],
       [-0.76341895, -0.21687009, -0.60840683],
       [-0.5186182 , -0.35564986,  0.77752707]])

The main problem I have is that the values for Principal Components are too different between each method. I would expect U and eigen_vectors to be close to each other. Am I missing something?

TLDR. Assuming that U, and eigen_vectors represent principal components for the original matrix X, why are they not equal?

Improve this question
$\endgroup$
11
  • 1
    $\begingroup$ It's not clear what you're specifically asking about. The eigenvalues are numerical. The linked thread discusses why the sign of the eigenvectors is arbitrary (if $v$ is an eigenvector to a matrix, then $-v$ is also). If you have a different question, please edit to explain what you know and what you would like to know. $\endgroup$
    – Sycorax ♦
    Commented Jan 14, 2022 at 20:36
  • $\begingroup$ @Sycorax Thanks, for the review. I've modified the question to make it easier to answer. Hope this helps now $\endgroup$
    – Alex.Kh
    Commented Jan 14, 2022 at 20:40
  • $\begingroup$ @Sycorax That's not quite true. Notice that the squares of the singular values of the column-centered matrix exactly equal (to numerical precision) twice ($=3-1$) the eigenvalues of the covariance matrix. This relationship holds generally for rectangular matrices. The issue in the question is that of course the eigenvalues of the square of a matrix (which is what the covariance is proportional to) will not equal the singular values of that matrix, as is evident for almost all $1\times 1$ matrices. The svd method in Python also returns the transpose of the eigenvectors. $\endgroup$
    – whuber ♦
    Commented Jan 14, 2022 at 20:51
  • 1
    $\begingroup$ Yes, I agree. I was reacting to an edit: OP wrote in a (removed) edit that they found $-V^T$ to match eigenvectors, which prompted my remark. If we take $X = USV^T$, one covariance eigen factorization is $\Sigma \propto X^TX=VS^2 V^T$. But if we expect the eigenvectors to be $U$, then something has gone wrong; perhaps we've transposed $X$. That's what I was trying to point out. (And likewise there's a similar discussion to be had about scaling factors...) This is why I suggest using a non-square $X$: you'll have the wrong number of vectors when you've transposed by mistake! $\endgroup$
    – Sycorax ♦
    Commented Jan 14, 2022 at 21:02
  • 1
    $\begingroup$ @Sycorax Thank you for the explanation. (I never saw that edit you were responding to...) $\endgroup$
    – whuber ♦
    Commented Jan 14, 2022 at 21:06

0

Reset to default