2

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?

amoeba
  • 93,463
  • 28
  • 275
  • 317
AliMM
  • 21
  • 1
  • I would guess it's because the scores that you get from `PcaHubert` are not 100% uncorrelated, as you can see by calling `cor(robtest@scores)`. I'd say that method 1 makes more sense to me in this case, because you get scores directly from the rotated loadings. – amoeba Jan 12 '18 at 07:47

0 Answers0