Eigenvalues of covariance matrix are negative I'm working on the PCA of the mnist dataset, and I get a very strange result, I created a matrix whose rows are flattened mnist images, When I try to compute the eigenvalues of the covariance matrix, I get some negative values. But the covariance matrix is positive semi-definite.
np.linalg.eigvals(np.dot(mnistBis[:, 0:20].T, mnistBis[:, 0:20])) # mnistBis.shape=(60000, 784)

array([ 4.79599869e+02, -1.19628465e+02,  9.68398702e+01,  1.88726171e-01,
        0.00000000e+00,  0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
        0.00000000e+00,  0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
        0.00000000e+00,  0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
        0.00000000e+00,  0.00000000e+00,  0.00000000e+00,  0.00000000e+00])

Modification:
Here is the full code
import tensorflow.keras.datasets.mnist as mnist
import matplotlib.pyplot as plt
import numpy as np

(x_train, y_train), (x_test, y_test) = mnist.load_data()
mnistBis = np.reshape(x_train, (-1, 28*28))
np.linalg.eigvalsh(np.dot(mnistBis[:, 0:20].T, mnistBis[:, 0:20]))

```
You seem to be using a generic eignvalue solver. Perhaps your matrix is badly conditioned and the aglorithm cannot pick it up. Lots of zero eigenvalues with few massive ones certainly seem suspicious. Have you tried solvers specific for hermitian/symmetric matricies? eigh &amp; eigvalsh
Also, I would z-score all your features, and drop all features with zero variance. This will set an upper bound on maximum possible eigenvalue since trace of your matrix will be the number of the features, and since positive definite matrix has only positive eigenvalues, the largest eigenvalue will be no larger than the trace.