# Sample R code to illustrate factor analysis
#
# Perform factor analysis on exams data
#
# step 1: load data 
# (as usual, adapt path to your personal configuration)
#
exams<-matrix(scan(file='exams.txt'),byrow=T,ncol=5)
nr<-length(exams[,1])
dimnames(exams)<-list(1:nr,c("vec","mec","alg","ala","sta"))
#
# step 2: factor analysis with one factor
# 
# Bartlett scores are eqn. (4.10) of online notes
# regression (Thomson) scores are eqn. (4.11) of online notes
#
exams.fa<-factanal(exams,factors=1,scores="regression")
print(exams.fa)
sco1<-exams.fa$scores
exams.fa<-factanal(exams,factors=1,scores="Bartlett")
sco2<-exams.fa$scores
# NB sco2=1.108076*sco1
plot(1:88,sco1)
#
# step 3: factor analysis with two factors
# 
exams.fa<-factanal(exams,factors=2,scores="regression")
print(exams.fa)
sco1<-exams.fa$scores
plot(1:88,sco1[,1])
plot(1:88,sco1[,2])
#
# produce biplot
#
biplot(sco1,exams.fa$loadings,cex=0.5)
#
# Factor analysis for tree ring data
#
T<-scan("2004GL021750-NOAMER.1400.txt",skip=1)
T<-matrix(T,nrow=581,ncol=72,byrow=TRUE)
tree<-T[1:581,3:72]
tt<- c(1400:1980)
#
tree1.fa<-factanal(tree,factors=1,scores="regression")
print(tree1.fa)
plot(tt,tree1.fa$scores)
#
tree2.fa<-factanal(tree,factors=2,scores="regression")
print(tree2.fa)
plot(tt,tree2.fa$scores[,1])
plot(tt,tree2.fa$scores[,2])
biplot(tree2.fa$scores,tree2.fa$loadings,cex=0.4)