01-08-Stat-Methods-Principle-Components.Rmd

## Principal Components {-}

```{r a24, comment=NA, warning=FALSE, message=FALSE, fig.width=8, fig.height=5}

## Generate some data
library(mvtnorm)

mu = rep(10, 6)
cov = matrix(nrow = 6, byrow = TRUE,
  c(1, 0, 1, 0, 1, 1,
    0, 2, 1, 2, 1, 1,
    1, 1, 2, 1, 1, 0,
    0, 2, 1, 3, 2, 1,
    1, 1, 1, 2, 4, 1,
    1, 1, 0, 1, 1, 10)
  )

set.seed(1000)
X = rmvnorm(100, mu, cov)

colnames(X) = paste("x", 1:6, sep = "")

## How do the data correlate?
round(cor(X), 3)
pairs(X)

summary(X)

## Build the principal components. We do not need to standardize the data
## since all of the variables are roughly the same scale.
(pr = prcomp(X))

## Principal components rotation matrix is actually the same as the eigen vectors
eigen(cov(X))$vectors

## Build the principal components from the coefficients
X.pca = X %*% pr$rotation

## The variance of the principal components is equal to the eigen values
eigen(cov(X))$values
diag(var(X.pca))

## PCA Summary
## The first two PCA account for 80% of all variation in the data
summary(pr)

X.pca = data.frame(X.pca)

## plot of the first two PCA
plot(X.pca$PC1, X.pca$PC2)

```