# How to perform dimensionality reduction with PCA in R

I have a big dataset and I want to perform a dimensionality reduction.

Now everywhere I read that I can use PCA for this. However, I still don't seem to get what to do after calculating/performing the PCA. In R this is easily done with the command princomp.

But what to do after calculating the PCA? If I decided I want to use the first $100$ principal components, how do I reduce my dataset exactly?

I believe what you are getting at in your question concerns data truncation using a smaller number of principal components (PC). For such operations, I think the function prcompis more illustrative in that it is easier to visualize the matrix multiplication used in reconstruction.

First, give a synthetic dataset, Xt, you perform the PCA (typically you would center samples in order to describe PC's relating to a covariance matrix:

#Generate data
m=50
n=100
frac.gaps <- 0.5 # the fraction of data with NaNs
N.S.ratio <- 0.25 # the Noise to Signal ratio for adding noise to data

x <- (seq(m)*2*pi)/m
t <- (seq(n)*2*pi)/n

#True field
Xt <-
outer(sin(x), sin(t)) +
outer(sin(2.1*x), sin(2.1*t)) +
outer(sin(3.1*x), sin(3.1*t)) +
outer(tanh(x), cos(t)) +
outer(tanh(2*x), cos(2.1*t)) +
outer(tanh(4*x), cos(0.1*t)) +
outer(tanh(2.4*x), cos(1.1*t)) +
tanh(outer(x, t, FUN="+")) +
tanh(outer(x, 2*t, FUN="+"))

Xt <- t(Xt)

#PCA
res <- prcomp(Xt, center = TRUE, scale = FALSE)
names(res)


In the results or prcomp, you can see the PC's (res$x), the eigenvalues (res$sdev) giving information on the magnitude of each PC, and the loadings (res$rotation). res$sdev
length(res$sdev) res$rotation
dim(res$rotation) res$x
dim(res$x)  By squaring the eigenvalues, you get the variance explained by each PC: plot(cumsum(res$sdev^2/sum(res$sdev^2))) #cumulative explained variance  Finally, you can create a truncated version of your data by using only the leading (important) PCs: pc.use <- 3 # explains 93% of variance trunc <- res$x[,1:pc.use] %*% t(res$rotation[,1:pc.use]) #and add the center (and re-scale) back to data if(res$scale != FALSE){
trunc <- scale(trunc, center = FALSE , scale=1/res$scale) } if(res$center != FALSE){
trunc <- scale(trunc, center = -1 * res$center, scale=FALSE) } dim(trunc); dim(Xt)  You can see that the result is a slightly smoother data matrix, with small scale features filtered out: RAN <- range(cbind(Xt, trunc)) BREAKS <- seq(RAN[1], RAN[2],,100) COLS <- rainbow(length(BREAKS)-1) par(mfcol=c(1,2), mar=c(1,1,2,1)) image(Xt, main="Original matrix", xlab="", ylab="", xaxt="n", yaxt="n", breaks=BREAKS, col=COLS) box() image(trunc, main="Truncated matrix (3 PCs)", xlab="", ylab="", xaxt="n", yaxt="n", breaks=BREAKS, col=COLS) box()  And here is a very basic approach that you can do outside of the prcomp function: #alternate approach Xt.cen <- scale(Xt, center=TRUE, scale=FALSE) C <- cov(Xt.cen, use="pair") E <- svd(C) A <- Xt.cen %*% E$u

#To remove units from principal components (A)
#function for the exponent of a matrix
"%^%" <- function(S, power)
with(eigen(S), vectors %*% (values^power * t(vectors)))
Asc <- A %*% (diag(E$d) %^% -0.5) # scaled principal components #Relationship between eigenvalues from both approaches plot(res$sdev^2, E$d) #PCA via a covariance matrix - the eigenvalues now hold variance, not stdev abline(0,1) # same results  Now, deciding which PCs to retain is a separate question - one that I was interested in a while back. Hope that helps. • Marc, you don't need to record center and scale explicitly, prcomp does that for you. Have a look at res$center and res$scale. IMHO it is less error prone to use these (no accidental difference about centering or not / scaling or not between the explicit call to scale and the prcomp call). Apr 28, 2013 at 9:58 • This answer needs expansion because it does not answer the OP's questions about what to do after calculating the PCA or how do I reduce my dataset exactly? Given that the OP has conducted PCA on his sample, his question is what to do with it, and what actually goes on with this subsamples; not how to do PCA. We might as well propose doing E <- eigen(cov(Sample)); A<- scale(scale=F, Sample) %*% E$vectors to get yet another way to get the scores also (which is actually what princomp does stats:::princomp.default). Apr 28, 2013 at 14:12
• @user11852 - the question specifically makes reference to reduction of the data set (i.e. the truncation that I have demonstrated here). I'll leave it up to him to decide whether or not this was what he was looking for. Apr 28, 2013 at 20:24
• @Marc, thanks for the response. I think I might need to step back and re-read everything again, because I am stuck on how any of the answer above deals with dimensionality reduction. Because as you show, dim(trunc) = dim(Xt). What was the benefit of it, the dimenions did not get reduced. Jul 18, 2013 at 18:04
• @B_Miner - Just keep in mind that truncation is used to focus on the main patterns in data and to filter out small scale patterns and noise. The truncated data is not smaller in terms of it's dimensions, but "cleaner". However, truncation does reduce the amount of data in that the entire matrix can be reconstructed with just a few vectors. A nice example is in the use of PCA for image compression, where a smaller number of PCs can be used to reconstruct the image. This smaller subset of vectors takes up less memory, but the reconstruction will have some loss in small scale detail. Jul 18, 2013 at 19:51

These other answers are very good and detailed, but I'm wondering if you're actually asking a much more basic question: what do you do once you have your PCs?

Each PC simply becomes a new variable. Say PC1 accounts for 60 % of the total variation and PC2 accounts for 30 %. As that's 90 % of the total variation accounted for, you could simply take these two new variables (PCs) as a simplified version of your original variables. This means fitting them to models, if that's what you're interested in. When it comes time to interpret your results, you do so in the context of the original variables that are correlated with each PC.

Sorry if I've underestimated the scope of the question!

I believe your original question stems from being a bit uncertain about what PCA is doing. Principal Component Analysis allows you to identify the principal mode of variation in your sample. Those modes are emperically calculated as the eigenvectors of your sample's covariance matrix (the "loadings"). Subsequently those vectors serve as the new "coordinate system" of your sample as you project your original sample in the space they define (the "scores"). The proportion of variation associated with the $i$-th eigenvector/ mode of variation/ loading/ principal component equals $\frac{\lambda_i}{\Sigma_{k=1}^{p} \lambda_k}$ where $p$ is your sample's original dimensionality ($p =784$ in your case). [Remember because your covariance matrix is non-negative Definite you'll have no negative eigenvalues $\lambda$.] Now, by definition the eigenvectors are orthogonal to each other. That means that their respective projections are also orthogonal and where originally you had a big possibly correlated sample of variables, now you have a (hopefully significantly) smaller linearly independent sample (the "scores").

Practically with PCA you are using the projections of the PCs (the "scores") as surrogate data for your original sample. You do all your analysis on the scores, and afterwards you reconstruct your original sample back using the PCs to find out out what happened on your original space (that's basically Principal Component Regression). Clearly, if you are able to meaningful interpreter your eigenvectors ("loadings") then you are in an even better position: You can describe what happens to your sample in the mode of variation presented by that loading by doing inference on that loading directly and not care about reconstruction at all. :)

In general what do you "after calculating the PCA" depends on the target of your analysis. PCA just gives you a linearly independent sub-sample of your data that is the optimal under an RSS reconstruction criterion. You might use it for classification, or regression, or both, or as I mentioned you might want to recognise meaningful orthogonal modes of variations in your sample.

A comment : I think the best naive way to decide the number of components to retain is to base your estimate on some threshold of sample variation you would like to retain in your reduced dimensionality sample rather than just some arbitrary number eg. 3, 100, 200. As user4959 explained you can check that cumulative variation by checking the relevant field of the list under the $loadingsfield in the list object produced by princomp. • As you mention Principal Component Regression, in R that is provided by package pls. As for the number of components to retain, I don't see any real advantage of deciding % variance over no. of components (maybe that is because I work with data that have very different levels of noise. As @Marc-in-the-box mentions, there are lots of different approaches to determine an appropriate no. of PCs, and the strategy will (should) depend both on the type data and on the type of data analysis that is to follow. Apr 28, 2013 at 10:04 • I said naive way; no scree plots, no average eigenvalue, no log-eigenvalue diagram, no partial correlation test. No probabilistic framework for your model. Because it is not what the OP asks. (I didn't mention the pls package by the way, princomp {stats}) There is a clear advantage of using % of cum. variance over a single arbitrary number: You know the quality of your reconstruction. If you just want to do dim. reduction to$K$dimensions ($K$<$D$,$D\$ being your original sample dimensionality) that's fine too but you don't really know the quality of your dimensionality reduction. Apr 28, 2013 at 13:28
• yes, it is naive. I don't agree that setting an arbitrary % of explained variance has an intrinsic advantage of any other arbitrary cut-off. But it's anyways not worth fighting about this as a) the OP never asked for advice about choosing the no. of PCs to keep and b) I think we agree that a proper inspection of the PCA model should be done anyways. Apr 28, 2013 at 15:37
• No problem; it was only a comment that I made before my answer anyway. (I'll put my comment paragraph last as I think it perplexes rather than clarifies what I want to say) Apr 28, 2013 at 15:44

After doing the PCA then you may select the first two components and plot.. You can see the variation of the components using a scree plot in R. Also using summary function with loadings=T you can fins the variation of features with the components.

Try to think what you want. You can interpret lots of things from PCA analysis.

Best Abhik