#=================================================================================#
#
# Model setup in UMX for multigroup.
#
#=================================================================================#


library(OpenMx)
require(umx)

# Read data
handcode <- c(1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 2, 1,
              1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 2, 2, 2, 2, 1, 
              1, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1, 1, 
              2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 
              2, 2, 2, 2, 1, 1, 1, 2, 1, 2, 1, 1, 2, 2)
alltask <- read.csv('SEM_data.csv')

mylabels<-names(alltask)

#=================================================================================#
# Data by group
RH_row <- which(handcode == 1)
dataRaw_RH <- mxData(observed = alltask[RH_row, ], type = 'raw')

LH_row <- which(handcode == 2)
dataRaw_LH <- mxData( observed=alltask[LH_row, ], type = "raw" )

#=================================================================================#
# Define Model1_RH
Model1RH <- umxRAM("Model1RH", data=dataRaw_RH, suffix = "_RH",
                   # Factor loadings
                   umxPath("F1", to = mylabels, labels = paste0("F1_to_",1:6), free = c(F,T,T,T,T,T), values = rep(1,6)),
                   umxPath("F2",to = mylabels, labels = paste0("F2_to_",1:6), free = c(F,T,T,T,T,T), values = c(0,rep(1,5))),
                   # Factor covariances: fixed to zero 
                   umxPath(cov = c("F1", "F2"),fixedAt=0),
                   # Factor variances: F1 free, F2 fixed at 1
                   umxPath(var='F1', freeAt=1),
                   umxPath(var='F2', fixedAt=1),
                   # Residual variances
                   umxPath(var=mylabels, labels = paste0("var_",1:6), freeAt=1),
                   # Means 
                   umxPath(means=mylabels, labels = paste0("mean_",1:6), freeAt=1),
                   umxPath(means=c('F1','F2'),fixedAt = 0))

Model1RH = mxRun(Model1RH)
umxSummary(Model1RH, std=T)
mxAPA(umxFitIndices(Model1RH))
plot(Model1RH,mean=FALSE)

#=================================================================================#

# Define Model1_LH
Model1LH <- umxRAM("Model1LH", data=dataRaw_LH, suffix = "_LH",
                   # Factor loadings
                   umxPath("F1", to = mylabels, labels = paste0("F1_to_",1:6), free = c(F,T,T,T,T,T), values = rep(1,6)),
                   umxPath("F2",to = mylabels, labels = paste0("F2_to_",1:6), free = c(F,T,T,T,T,T), values = c(0,rep(1,5))),
                   # Factor covariances: fixed to zero 
                   umxPath(cov = c("F1", "F2"),fixedAt=0),
                   # Factor variances: F1 free, F2 fixed at 1
                   umxPath(var='F1', freeAt=1),
                   umxPath(var='F2', fixedAt=1),
                   # Residual variances
                   umxPath(var=mylabels, labels = paste0("var_",1:6), freeAt=1),
                   # Means 
                   umxPath(means=mylabels, labels = paste0("mean_",1:6), freeAt=1),
                   umxPath(means=c('F1','F2'),fixedAt = 0))

Model1LH = mxRun(Model1LH)
umxSummary(Model1LH, std=T)
umxAPA(umxFitIndices(Model1LH))
plot(Model1LH,mean=FALSE)


#=================================================================================#

# Multigroup modelling using UMX

mgFitFun <- umxSuperModel("myfirstsupermodel", Model1RH, Model1LH)

umxSummary(mgFitFun,std=T)

#=================================================================================#