Let me continue my comment with an illustration for the case where you're interested in existing R packages. There are several package in the Multivariate Task View that will provide enhanced method for PCA-related methods (as compared to R base prcomp
and princomp
), e.g. ade4 or FactorMineR. I personally like FactoMineR because of its simple syntax, and you check the associated website for more information on the available methods.
One can use supplementary categorical and/or numerical variables when applying a PCA. Those variables are not used to construct factor axes, but can be showed afterwards on the correlation circle (for numerical variables) or the individual map (for categorical variables). Here is a toy example of use (from the on-line help):
data(decathlon)
res.pca <- PCA(decathlon, quanti.sup = 11:12, quali.sup=13)
plotellipses(res.pca,13)
If you have multiple passive variables, you can select the one to display (with or without confidence ellipses) using the keepvar=
argument. Here is another picture with two illustrative variables.
Be careful with arguments that are a little bit non-standard if you are used to default plotting functions in R. The plotellipses()
function makes use of the helper function ellipse::ellipse
that you can use (or not) in any plot (look for monpanel.ellipse
subfunction in plotellipses()
to see how confidence lines are computed). That's what I did to build specific individual map (B&W, different plotting symbol, etc.). For example, the following snippet just plot all individuals with two different symbols depending on the type of sporting event (2004 Olympic Game or 2004 Decastar):
labs <- paste(round(res.pca$eig[1:2, 2], 2), "%", sep="")
plot(res.pca$ind$coord[,1:2], pch=as.numeric(decathlon$Competition),
xlab=paste("Dim. 1 (", labs[1], ")", sep=""),
ylab=paste("Dim. 2 (", labs[2], ")", sep=""))
abline(v=0, h=0, lty=2)
Besides, I would like to point you to @vqv's excellent ggbiplot package, available on GitHub, which follows from one of his answer. (It uses R base functions and ggplot2.)
mydata
, hence your workingeset
. Any chance we can get simulated data or a pointer to an existing R/bioc dataset? $\endgroup$