library(haven)
SCoA_HigherEd_WLSMV <- read_sav("C:/Users/gbro008/Google Drive/2019_08 Moscow summer school/Lectures/resources/SCoA HigherEd WLSMV.sav")
View(SCoA_HigherEd_WLSMV)


#install required libraries
library(lavaan)
library(semPlot)
#to create 8 factor SCoA model
SCoA6.modelwls <- '
Bad =~ bd1 + bd2 + bd3 + bd4 + bd5
CE =~ ce1 + ce2 + ce3 + ce4 + ce5 + ce6
IG =~ ig1 + ig2 + ig3
PE =~ pe1 + pe2
SF =~ sf1 + sf2 + sf3
SI =~ si1 + si2 + si3 + si4 + si5
SQ =~ sq1 + sq2
TI =~ ti1 + ti2 + ti3 + ti4 + ti5 + ti6
'


#cfa analysis with WLSmv estimator, ordered command to inform CFA that this is ordered scores
SCoAwls_fit <-cfa (SCoA6.modelwls, data=SCoa_HErand, ordered=c("bd1", "bd2", "bd3", "bd4", "bd5", "ce1", "ce2", "ce3", "ce4", "ce5", "ce6", "ig1", "ig2", "ig3", "pe1", "pe2", "sf1", "sf2", "sf3", "si1", "si2", "si3", "si4", "si5", "sq1", "sq2", "ti1", "ti2", "ti3", "ti4", "ti5", "ti6"))



#get results of cfa
summary (SCoAwls_fit, standardized=TRUE)
#detailed fit indices
fitmeasures(SCoAwls_fit)
#to get gamma hat fit index
library(semTools)
moreFitIndices(SCoAwls_fit)
#to get x2/df ratio and p-value
fm <- fitMeasures(SCoAwls_fit)
fm[["chisq"]] / fm[["df"]]
fm["pvalue"] 
#to get conventional correlated factors diagram
semPaths(SCoAwls_fit, intercept = FALSE, whatLabel = "std",  residuals = FALSE, exoCov = TRUE)