Computing SVD through Eigendecomposition of correlation matrix I am following the excellent series on SVD by Steve Brunton from the University of Washington, on YouTube, but I have trouble interpreting his 4th video on the subject.
If I understand correctly, he mentions that one can compute the economy SVD decomposition $X = \hat{U}\hat{\Sigma}V^T$ with the following :
$$X^TX= V\hat{\Sigma}\hat{U}^T\hat{U}\hat{\Sigma}V^T = V\hat{\Sigma}^2V^T \implies X^TXV = V\hat{\Sigma}^2 $$
$$XX^T= \hat{U}\hat{\Sigma}\hat{V}^T\hat{V}\hat{\Sigma}\hat{U}^T = \hat{U}\hat{\Sigma}^2\hat{U}^T \implies XX^T\hat{U} = \hat{U}\hat{\Sigma}^2 $$
Instead of using U,S,VT = svd(X) in Python, I want to try decomposing the image in its SVD, and reconstructing the original using only $r$ values with this technique. I am trying to apply this on a picture of The Starry Night, loaded in a grayscale numpy array X, which I do this way :
r = 10    
XT = X.transpose()
C_1 = XT @ X
C_2 = X @ XT
[Lambda_1, V_hat] = np.linalg.eig(C_1)
[Lambda_2, U_hat] = np.linalg.eig(C_2)
V_hat_T = V_hat.transpose()
X_tilde = U_hat[:,:r] @ np.diag(Lambda_1)[0:r,:r] @ V_hat_T[:r,:]

But I can't reproduce the image, whatever the value I attribute to $r$, whereas the product of U @ S @ VT , computed with the svd(X) function, perfectly reconstructs the original picture.
What am I doing wrong? It is certainly superfluous to mention that I am a beginner, and I probably made some big mistake.
The problem is occurring for a few reasons:

You must ensure the correct ordering of eigenvalues/eigenvectors.
You should be taking the square root of the eigenvalues (Lambda_1) to get the singular values.
You need to take special care to ensure you are using the right signs for the eigenvectors as these won't necessarily be correct by default. This is because we assume that the singular values are the positive roots of the eigenvalues, but the two sets of eigenvectors we retrieve don't necessarily respect this assumption when joined in the singular value decomposition.

Here is the working example:
r = 10    
XT = X.transpose()
C_1 = XT @ X
C_2 = X @ XT
[Lambda_1, V_hat] = np.linalg.eig(C_1)
[Lambda_2, U_hat] = np.linalg.eig(C_2)

# The eigenvalues/eigenvectors returned from np.linalg.eig are not sorted so may differ between Lambda_1 and Lambda_2
# When performing economy SVD we also want to ensure the largest eigenvalues come first
i_1 = np.argsort(Lambda_1)[::-1]
i_2 = np.argsort(Lambda_2)[::-1]
V_hat = V_hat[:, i_1]
U_hat = U_hat[:, i_2]

# We must take the square root of eigenvalues to get the singular values
Lambda = np.sqrt(np.sort(Lambda_1)[::-1])

# X @ V_hat and U_hat @ np.diag(Lambda) should be equal. But the eigenvectors retrieved may not have the correct signs 
# This code checks if  X @ V_hat and U_hat @ np.diag(Lambda) have different signs and then updates V_hat wherever the sign is incorrect
same_sign = np.sign((X @ V_hat)[0] * (U_hat @ np.diag(Lambda))[0])
V_hat = V_hat * same_sign.reshape(1, -1)

V_hat_T = V_hat.transpose()
X_tilde = U_hat[:,:r] @ np.diag(Lambda)[0:r,:r] @ V_hat_T[:r,:]