This great answer shows how to compute varimax-rotated loadings and scores from PCA results using princomp.
I am conducting a robust PCA analysis using the pcaHubert function from the rrcov package. While I can rotate the loadings without any issues, when I try to reproduce the rotated scores using the second and third methods from the linked answer above, I get slightly different results between the two methods. I am not sure why, but guess it must be something related to differences the robust/ROBPCA output or algorithm? Can someone explain what is causing this, how to fix it, and/or if one of these methods is "more' correct in the case of robust PCA?
irisX <- iris[,1:4] #get example data ncomp <- 2 #number of PCs robtest <- PcaHubert(irisX, k=ncomp, alpha =0.75, scale=TRUE) #Robust pca #Rotate loadings sdev <- sqrt(robtest@eigenvalues) #compute SD of PCs from eigenvalues rawload <- robtest@loadings[,1:2] %*% diag(sdev, ncomp, ncomp) #convert to raw loadings vmxload <- varimax(rawload)$loadings #varimax rotated loadings #Rotate scores #Method 1 invLoadings <- t(pracma::pinv(vmxload)) scores1 <- scale(irisX) %*% invLoadings #multiply data with the transposed pseudo-inverse of the rotated loadings #Method2 scores2 <- scale(robtest@scores[,1:2]) %*% varimax(rawload)$rotmat #Compare results > print(scores1[1:2,]) [,1] [,2] [1,] 1.024056 -0.8049907 [2,] 1.365848 0.3179221 > print(scores2[1:2,]) [,1] [,2] [1,] 1.093374 -0.8408940 [2,] 1.433536 0.3005973
Clearly I'm missing something, maybe quite basic, but can't fathom what is is?