#GO TO http://cran.r-project.org/ AND DOWNLOAD THE R RELEASE: R-2.5.1.
# or http://www.staff.ul.ie/mackenzieg/Assess/R_Integration_Tools/r_integration_tools.html
# for a 'crash course' in R

#PUT THE SPSS FILE ON A MEMORY STICK IN THE BIN FOLDER OF R-2.5.1.
#MAKE SURE YOU HAVE LIBRARY FOLDERS FOREIGN, MVTNORM, SFSMISC AND POLYCOR IN THE
#R-2.5.1 LIBRARY FOLDER; THESE MAY BE DOWNLOADED FROM http://cran.r-project.org/
#AND CLICKING ON CONTRIBUTED EXTENSION PACKAGES.

#(PW: USED R-2.5.1 ON BLACK MEMORY STICK; ONLY WORKS WITH R-2.5)
#( CAN THEN COPY CORRELATIONS INTO A .DAT FILE AND READ INTO SPSS; CREATE SEPARATE FILES AND COPY AND PASTE;
#THEN ENTER INTO A FACTOR ANALYSIS)

#(AMEND SPSS FILENAME AND ITEM NAMES IN THE BELOW R CODE WHICH ARE ASSUMED AS SATCAT, DAT, PERF_2 AND PERF_3)

#CLICK RGUI.EXE AND CUT AND PASTE IN THE AMENDED CODE:
#
library(foreign)
library(mvtnorm)
library(sfsmisc)
library(polycor)
x <- read.spss("E:\\FACTOR FOR ZIPPING\\DEMO FACTOR\\example33.sav")
p1 <- as.factor(x$ptci1)
p2 <- as.factor(x$ptci2)
p3 <- as.factor(x$ptci3)
p4 <- as.factor(x$ptci4)
p5 <- as.factor(x$ptci5)
p6 <- as.factor(x$ptci6)
res <- hetcor(p1, p2, p3, p4, p5, p6, ML=FALSE, std.err=FALSE)
res$correlation
library(psych)
KMO(res$correlation)
res$type
res$correlation[upper.tri(res$correlation)] <- 0
res1 <- factanal(covmat=res$correlation, factors=3, rotation="varimax")
res2 <- factanal(covmat=res$correlation, factors=3, rotation="promax")
res2 <- factanal(covmat=res$correlation, factors=3, rotation="none")

res3 <- princomp(covmat=res$correlation)
pdf("scree_princomp.pdf")
screeplot(res3)
dev.off()
library(psych)
xdat <- cbind(p1,p2,p3,p4,p5,p6)
xdat <- data.frame(xdat)
res3 <- principal(xdat,rotate="varimax",nfactors=3)
pdf("E:\\FACTOR FOR ZIPPING\\DEMO FACTOR\\scree_principal.pdf")
plot(res3$values,type="b")
dev.off()

# PW NOTE: SOME MISSING LOADINGS; THE MISSINGNESS PATTERN DEPENDS ON WHICH OF THE TWO ROTATIONS IS USED;

library(sem)
model.f1<-matrix(c(    'F1->p1', 'lam1',NA,
                       'F1->p2', 'lam2',NA,
                       'F1->p3', 'lam3',NA,
                       'F1->p4', 'lam4',NA,
                       'F1->p5', 'lam5',NA,
                       'F1->p6', 'lam6',NA,
                      'p1<->p1', 'th1' ,NA,
                      'p2<->p2', 'th2' ,NA,
                      'p3<->p3', 'th3' ,NA,
                      'p4<->p4', 'th4' ,NA,
                      'p5<->p5', 'th5' ,NA,
                      'p6<->p6', 'th6' ,NA,
                      'F1<->F1',   NA , 1),
             ncol=3, byrow=T)
sem.poly.f1<-sem(model.f1, res$correlation, 278)
summary(sem.poly.f1)

model.f1<-matrix(c(    'F1->p1', 'lam1',NA,
                       'F1->p2', 'lam2',NA,
                       'F1->p3', 'lam3',NA,
                       'F2->p4', 'lam4',NA,
                       'F2->p5', 'lam5',NA,
                       'F2->p6', 'lam6',NA,
                      'p1<->p1', 'th1' ,NA,
                      'p2<->p2', 'th2' ,NA,
                      'p3<->p3', 'th3' ,NA,
                      'p4<->p4', 'th4' ,NA,
                      'p5<->p5', 'th5' ,NA,
                      'p6<->p6', 'th6' ,NA,
                      'F1<->F1',   NA,  1,
                      'F2<->F2',   NA,  1,
                      'F2<->F1', 'gam1', NA),
             ncol=3, byrow=T)
sem.poly.f1<-sem(model.f1, res$correlation, 278)
summary(sem.poly.f1)

model.f1<-matrix(c(    'F1->p1', 'lam1',NA,
                       'F1->p2', 'lam2',NA,
                       'F1->p3', 'lam3',NA,
                       'F1->p5', 'lam5',NA,
                       'F2->p1', 'lam6',NA,
                       'F2->p2', 'lam7',NA,
                       'F2->p4', 'lam9',NA,
                       'F2->p6', 'lam10',NA,
                      'p1<->p1', 'th1' ,NA,
                      'p2<->p2', 'th2' ,NA,
                      'p3<->p3', 'th3' ,NA,
                      'p4<->p4', 'th4' ,NA,
                      'p5<->p5', 'th5' ,NA,
                      'p6<->p6', 'th6' ,NA,
                      'F1<->F1',   NA,  1,
                      'F2<->F2',   NA,  1),
             ncol=3, byrow=T)
sem.poly.f1<-sem(model.f1, res$correlation, 278)
summary(sem.poly.f1)

# ALGORITHM DOES NOT ALWAYS CONVERGE! COMMON PROBLEM IN SEM - NEED GOOD STARTING VALUES!

model.f1<-matrix(c(    'F1->p1', 'lam1',NA,
                       'F1->p2', 'lam2',NA,
                       'F1->p3', 'lam3',NA,
                       'F1->p4', 'lam5',NA,
                       'F1->p5', 'lam6',NA,
                       'F2->p1', 'lam7',NA,
                       'F2->p2', 'lam8',NA,
                       'F2->p3', 'lam9',NA,
                       'F2->p4', 'lam10',NA,
                       'F2->p6', 'lam11',NA,
                      'p1<->p1', 'th1' ,NA,
                      'p2<->p2', 'th2' ,NA,
                      'p3<->p3', 'th3' ,NA,
                      'p4<->p4', 'th4' ,NA,
                      'p5<->p5', 'th5' ,NA,
                      'p6<->p6', 'th6' ,NA,
                      'F1<->F1',   NA,  1,
                      'F2<->F2',   NA,  1),
             ncol=3, byrow=T)
sem.poly.f1<-sem(model.f1, res$correlation, 278)
summary(sem.poly.f1)

sem.poly.f1

# CAN USE DIFFERENT STARTING VALUES - IF YOU CAN FIND ONES THAT AID CONVERGENCE
# OR CHANGE MODEL - MAY BE A BAD FIT!

model.f1<-matrix(c(    'F1->p1', 'lam1',0.012,
                       'F1->p2', 'lam2',NA,
                       'F1->p3', 'lam3',NA,
                       'F1->p4', 'lam5',NA,
                       'F1->p5', 'lam6',NA,
                       'F2->p1', 'lam7',0.381,
                       'F2->p2', 'lam8',0.629,
                       'F2->p3', 'lam9',0.705,
                       'F2->p4', 'lam10',0.562,
                       'F2->p6', 'lam11',0.342,
                      'p1<->p1', 'th1' ,NA,
                      'p2<->p2', 'th2' ,NA,
                      'p3<->p3', 'th3' ,NA,
                      'p4<->p4', 'th4' ,NA,
                      'p5<->p5', 'th5' ,NA,
                      'p6<->p6', 'th6' ,NA,
                      'F1<->F1',   NA,  1,
                      'F2<->F2',   NA,  1),
             ncol=3, byrow=T)
sem.poly.f1<-sem(model.f1, res$correlation, 278)
summary(sem.poly.f1)