Does the new OpenMX version allow me to model a Four-factor Cholesky Decomposition? I would like to use 2 continuous, 1 ordinal, and 1 dichotomous variable in the model? I was wondering if anyone might have code that is similar to post on the forum for 3 or 4 factor models? Below is my attempt but I get NaN for the standard error.
vars v1 and v2 are continuous and v3 is ordinal (5 categories), but v4 is dichotomous
nvar <- 4 # number of variables
tnvar <-8 # number of variablesmax family size
nlower <-tnvar(tnvar+1)/2 # number of free elements in a lower matrix tnvar*tnvar
Vars <- c("v1","v2","v3","v4")
selVars <- paste(Vars,c(rep(1,nvar),rep(2,nvar)),sep="")
This yields a list of the 6 variables for analysis
v11 v21 v31 v41 v12 v22 v32 v42
mzData <- as.matrix(subset(octotwin_data, zygosity==1, selVars))
dzData <- as.matrix(subset(octotwin_data, zygosity==2, selVars))
Print Descriptive Statistics
-------------------------------------------------------------------------
summary(mzData)
summary(dzData)
Fit Multivariate ACE Model with RawData and Matrices Input
ACE_Cholesky_Model <- mxModel("ACE_Cholesky",
mxModel("ACE",
# Matrices a, c, and e to store a, c, and e path coefficients
mxMatrix( type="Lower", nrow=nvar, ncol=nvar, free=TRUE, values=.6, name="a" ),
mxMatrix( type="Lower", nrow=nvar, ncol=nvar, free=TRUE, values=.6, name="c" ),
mxMatrix( type="Lower", nrow=nvar, ncol=nvar, free=TRUE, values=.6, name="e" ),
# Matrices A, C, and E compute variance components
mxAlgebra( a %% t(a), name="A" ),
mxAlgebra( c %% t(c), name="C" ),
mxAlgebra( e %% t(e), name="E" ),
# Algebra to compute total variances and standard deviations (diagonal only)
mxAlgebra( A+C+E, name="V" ),
mxMatrix( type="Iden", nrow=nvar, ncol=nvar, name="I"),
mxAlgebra( solve(sqrt(IV)), name="iSD"),
# Matrix & Algebra for expected means vector
mxMatrix( type="Full", nrow=1, ncol=nvar, free=TRUE, values= 0, name="Mean" ),
mxAlgebra( cbind(Mean,Mean), name="expMean"),
# Algebra for expected variance/covariance matrix in MZ
mxAlgebra( rbind ( cbind(A+C+E , A+C),
cbind(A+C , A+C+E)), name="expCovMZ" ),
# Algebra for expected variance/covariance matrix in DZ
mxAlgebra( rbind( cbind(A+C+E , 0.5%x%A+C),
cbind(0.5%x%A+C , A+C+E)), name="expCovDZ" )
),
mxModel("MZ",
mxData( observed=mzData, type="raw" ),
mxFIMLObjective( covariance="ACE.expCovMZ", means="ACE.expMean", dimnames=selVars )
),
mxModel("DZ",
mxData( observed=dzData, type="raw" ),
mxFIMLObjective( covariance="ACE.expCovDZ", means="ACE.expMean", dimnames=selVars )
),
mxAlgebra( MZ.objective + DZ.objective, name="modelfit" ),
mxAlgebraObjective("modelfit")
)
ACE_Cholesky_Fit <- mxRun(ACE_Cholesky_Model)
ACE_Cholesky_Summ <- summary(ACE_Cholesky_Fit)
ACE_Cholesky_Summ
Generate Multivariate Cholesky ACE Output
parameterSpecifications(ACE_Cholesky_Fit)
expectedMeansCovariances(ACE_Cholesky_Fit)
tableFitStatistics(ACE_Cholesky_Fit)
ACEpathMatrices <- c("ACE.a","ACE.c","ACE.e","ACE.iSD","ACE.iSD %% ACE.a","ACE.iSD %% ACE.c","ACE.iSD %*% ACE.e")
ACEpathLabels <- c("pathEst_a","pathEst_c","pathEst_e","sd","stPathEst_a","stPathEst_c","stPathEst_e")
formatOutputMatrices(ACE_Cholesky_Fit,ACEpathMatrices,ACEpathLabels,Vars,3)
ACEcovMatrices <- c("ACE.A","ACE.C","ACE.E","ACE.V","ACE.A/ACE.V","ACE.C/ACE.V","ACE.E/ACE.V")
ACEcovLabels <- c("covComp_A","covComp_C","covComp_E","Var","stCovComp_A","stCovComp_C","stCovComp_E")
formatOutputMatrices(ACE_Cholesky_Fit,ACEcovMatrices,ACEcovLabels,Vars,3)
Genetic and Environmental Correlations
ACEcorMatrices <- c("solve(sqrt(ACE.IACE.A)) %% ACE.A %% solve(sqrt(ACE.IACE.A))",
"solve(sqrt(ACE.IACE.C)) %% ACE.C %% solve(sqrt(ACE.IACE.C))",
"solve(sqrt(ACE.IACE.E)) %% ACE.E %% solve(sqrt(ACE.IACE.E))",
"solve(sqrt(ACE.IACE.V)) %% ACE.V %% solve(sqrt(ACE.IACE.V))")
ACEcorLabels <- c("corA","corC","corE","cor")
formatOutputMatrices(ACE_Cholesky_Fit,ACEcorMatrices,ACEcorLabels, Vars, 4)