EU <- read.table("clipboard")

# 'EU' data set with 3 columns in data frame EU

par(mfrow=c(1,2))
plot(EU[,1], EU[,2], type="n", cex.axis=0.7, xlab="Purchasing power/capita", ylab="GDP/capita")
text(EU[,1], EU[,2], labels=rownames(EU), col="blue", font=2)
plot(EU[,1], EU[,3], type="n", cex.axis=0.7, xlab="Purchasing power/capita", ylab="Inflation rate")
text(EU[,1], EU[,3], labels=rownames(EU), col="blue", font=2)

# first install the R package rgl before carrying on

library(rgl)

plot3d(EU[,1], EU[,2], EU[,3], xlab="Purchasing power", ylab="GDP", zlab="Inflation", font=2, col="red", type="n")
text3d(EU[,1], EU[,2], EU[,3], text=rownames(EU), font=2, col="blue") 

for(ir in 1:181) {
  theta <- 360*(ir-1)/181
  view3d( theta = theta, phi = 15 , fov = 20, zoom=0.8)
  Sys.sleep(0.1) 
}

# correlations between the three variables

cor(EU)

# readership data
readers <- read.table("clipboard")
row.pro <- readers/apply(readers, 1, sum)

plot3d(row.pro[,1], row.pro[,3], row.pro[,2], xlab="C1", ylab="C3", zlab="C2", font=2, col="red", type="n")
text3d(row.pro[,1], row.pro[,3], row.pro[,2], text=rownames(row.pro), font=2, col="blue") 

row.pro <- rbind(as.matrix(row.pro), diag(c(1,1,1)))
rownames(row.pro)[6:8] <- colnames(row.pro)

plot3d(row.pro[,1], row.pro[,3], row.pro[,2], xlab="C1", ylab="C3", zlab="C2", font=2, col="red", type="n")
text3d(row.pro[,1], row.pro[,3], row.pro[,2], text=rownames(row.pro), font=2, col=c(rep("blue",5),rep("red",3))) 
lines3d(row.pro[c(6:8,6),1],row.pro[c(6:8,6),3],row.pro[c(6:8,6),2], col="red") 


# read in data into data-frame data_set
# the next 14 commands are all you need to compute CA results
data.P    <- data_set/sum(data_set)
data.r    <- apply(data.P,1,sum)
data.c    <- apply(data.P,2,sum)
data.Drmh <- diag(1/sqrt(data.r))
data.Dcmh <- diag(1/sqrt(data.c))
data.P    <- as.matrix(data.P)
data.S    <- data.Drmh %*% (data.P-data.r%o%data.c) %*% data.Dcmh
data.svd  <- svd(data.S)
data.rsc  <- data.Drmh%*%data.svd$u
data.csc  <- data.Dcmh%*%data.svd$v
data.rpc  <- data.rsc%*%diag(data.svd$d)
data.cpc  <- data.csc%*%diag(data.svd$d)
# the symmetric map
plot(c(data.rpc[,1],data.cpc[,1]),c(data.rpc[,2],data.cpc[,2]), type="n", xlab="dim1", ylab="dim2")
text(data.rpc[,1],data.rpc[,2],label=rownames(data_set),col="blue", font=2)
text(data.cpc[,1],data.cpc[,2],label=colnames(data_set), col="red", font=4)
# now do it in one shot using ca package (first install from CRAN)
library(ca)
plot(ca(data_set))


# compute matrix of contributions for rows and inertias
data.rcon <- data.rpc^2 * data.r
apply(data.rcon, 1, sum)

# compute contributions and squared correlations
data.rctr <- t( t(data.rcon) / apply(data.rcon, 2, sum) )
data.rcor <- data.rcon / apply(data.rcon, 1, sum)

# compute qualities in 2-d solution
apply(data.rcor[,1:2], 1, sum)

# compute matrix of contributions for columns and inertias
data.ccon <- data.cpc^2 * data.c
apply(data.ccon, 1, sum)

# compute contributions and squared correlations
data.cctr <- t( t(data.ccon) / apply(data.ccon, 2, sum) )
data.ccor <- data.ccon / apply(data.ccon, 1, sum)

# compute qualities in 2-d solution
apply(data.ccor[,1:2], 1, sum)


# biplot

# PCA of EU data (with standardization of variables)
EU.PCA<-princomp(EU, cor=T)
biplot(EU.PCA, scale=0)

# ca biplot -- asymmetric map and contribution biplot

data(author)

plot(ca(author), map="rowprincipal")
plot(ca(author), map="rowgreen")

# option for symbol to be proportional to mass

plot(ca(author), map="rowgreen", mass=c(F,T))

# 3d, with lines out to the column points

plot3d(ca(author), map="rowgreen", arrows=c(F,T), sf=0.000001)