# NOT RUN {
# Example data from ?prcomp (see discussion at Stats.StackExchange.com/q/397793)
C <- chol(S <- toeplitz(.9 ^ (0:31)))
set.seed(17)
X <- matrix(rnorm(32000), 1000, 32)
Z <- X %*% C
pcaz <- prcomp(Z)
tst <- PCAtoXhat(pcaz)
all.equal(tst, Z, check.attributes = FALSE)
# Plot to show the effect of increasing ncomp
ntests <- ncol(Z)
rmsd <- rep(NA_real_, ntests)
for (i in 1:ntests) {
ans <- XtoPCAtoXhat(X, i, sd)
del<- ans - X
rmsd[i] <- sqrt(sum(del^2)/length(del)) # RMSD
}
plot(rmsd, type = "b",
main = "Root Mean Squared Deviation\nReconstructed - Original Data",
xlab = "No. of Components Retained", ylab = "RMSD")
abline(h = 0.0, col = "pink")
# }
Run the code above in your browser using DataLab