CODEWITHSUNDEEP
rpca
Federico Giorgi
Federico Giorgi

Reputation: 10755

R: backwards principal component calculation

I would like to perform a backwards principal component calculation in R, meaning: obtaining the original matrix by the PCA object itself.

This is an example case:

# Load an expression matrix
load(url("http://www.giorgilab.org/allexp_rsn.rda"))
# Calculate PCA
pca <- prcomp(t(allexp_rsn))

In order to obtain the original matrix, one should multiply the rotations by the PCA themselves, as such:

test<-pca$rotation%*%pca$x

However, as you may check, the calculated "test" matrix is completely different from the original "allexp_rsn" matrix. What am I doing wrong? Is the function prcomp adding something else to the svs procedure?

Thanks :-)

Upvotes: 2

Views: 1048

Answers (3)

Thrawn
Thrawn

Reputation: 87

Here is a solution using the eigen function, applied to a B/W image matrix to illustrate the point. The function uses increasing numbers of PCs, but you can use all of them, or only some of them

library(gplots)
library(png) 
# Download an image:
download.file("http://www.giorgilab.org/pictures/monalisa.tar.gz",destfile="monalisa.tar.gz",cacheOK = FALSE)
untar("monalisa.tar.gz")
# Read image:
img <- readPNG("monalisa.png")
# Dimension
d<-1
# Rotate it:
rotate <- function(x) t(apply(x, 2, rev))
centermat<-rotate(img[,,d])
# Plot it
image(centermat,col=gray(c(0:100)/100))
# Increasing PCA
png("increasingPCA.png",width=2000,height=2000,pointsize=20)
par(mfrow=c(5,5),mar=c(0,0,0,0))
for(end in (1:25)*12){
    for(d in 1){
        centermat<-rotate(img[,,d])
        eig <- eigen(cov(centermat))
        n <- 1:end
        eigmat<-t(eig$vectors[,n] %*% (t(eig$vectors[,n]) %*% t(centermat)))
        image(eigmat,col=gray(c(0:100)/100))
    }
}
dev.off()

backwards PCA reconstruction by Federico Manuel Giorgi

Upvotes: 1

Banach
Banach

Reputation: 158

Remember that a prerequisite to perform PC analysis is to scale and center the data. I believe that prcomp procedure does that, so pca$x returns scaled original data (with mean 0 and std. equal to 1).

Upvotes: 1

user3710546
user3710546

Reputation:

Using USArrests:

pca <- prcomp(t(USArrests))
out <- t(pca$x%*%t(pca$rotation))
out <- sweep(out, 1, pca$center, '+')
apply(USArrests - out, 2, sum)
   Murder       Assault      UrbanPop          Rape 
1.070921e-12 -2.778222e-12  3.801404e-13  1.428191e-12

Upvotes: 2

Related Questions