#reading data corMZ <- read.csv("MZcor.csv") corDZ <- read.csv("DZcor.csv") lowerTriangle(corMZ1) <- lowerTriangle(t(corMZ1)) lowerTriangle(corDZ1) <- lowerTriangle(t(corDZ1)) #defining variables variables <- c("FAM1","FAM2","A1","A2","A3","A4","A5","B1","B2","B3","B4","B5") varsFAM <- c("FAM1","FAM2") varsA <- c("A1","A2","A3","A4","A5") varsB <- c("B1","B2","B3","B4","B5") # defining labels for openmx within <- paste("a",1:15,sep="") withinF <- paste("b",1:10,sep="") mz <- paste("mz",1:15,sep="") dz <- paste("dz",1:15,sep="") #you may need to run a line similar to the one below corMZ <- matrix(corMZ, nrow=12,dimnames = list(variables,variables)) corDZ <- matrix(corDZ, nrow=12,dimnames = list(variables,variables)) #rather insane data driven start values for openmx, there have been ~20 different tries on the start values withinS <- (upperTriangle(t(corDZ[varsA,varsA]),diag=T)+upperTriangle(t(corDZ[varsB,varsB]),diag=T)+upperTriangle(t(corMZ[varsA,varsA]),diag=T) +upperTriangle(t(corMZ[varsB,varsB]),diag=T))/4 withinFS <- c((corMZ[varsFAM,varsA]+corMZ[varsFAM,varsB]+corDZ[varsFAM,varsA]+corDZ[varsFAM,varsB])/4) acrossmS <- upperTriangle(t(corMZ[varsA,varsB]),diag=T) acrossdS <- upperTriangle(t(corDZ[varsA,varsB]),diag=T) withinFAMS <- upperTriangle(t(corDZ[varsFAM,varsFAM]),diag=T) #I only had one free parameter in the one variable that ends up being exogenous in our path diagram #since it's a correlation matrix, not sure if I should set all diag to fixed... ModelSAT <- mxModel("SAT", mxModel("MZ", type="RAM", manifestVars=vars, mxPath(from=varsA, to=varsA, arrows=2, free=T,connect="unique.pairs", values=withinS,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsB, to=varsB, arrows=2, free=T,connect="unique.pairs", values=withinS,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsA, arrows=2, free=T,connect="unique.pairs", values=withinFS,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsB, arrows=2, free=T,connect="unique.pairs", values=withinFS,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsFAM, arrows=2, free=c(F,T,T),connect="unique.pairs", values=withinFAMS,labels=c("c1","c2","c3"), lbound=-1,ubound=1.5), mxPath(from=varsA, to=varsB, arrows=2, free=T,connect="unique.pairs", values=acrossmS,labels=mz, lbound=-1,ubound=1.5), mxData( observed=corMZ, type="cor" ,numObs=203)), mxModel("DZ", type="RAM", manifestVars=vars, mxPath(from=varsA, to=varsA, arrows=2, free=T,connect="unique.pairs", values=withinS,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsB, to=varsB, arrows=2, free=T,connect="unique.pairs", values=withinS,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsA, arrows=2, free=T,connect="unique.pairs", values=withinFS,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsB, arrows=2, free=T,connect="unique.pairs", values=withinFS,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsFAM, arrows=2, free=c(F,T,T),connect="unique.pairs", values=withinFAMS,labels=c("c1","c2","c3"), lbound=-1,ubound=1.5), mxPath(from=varsA, to=varsB, arrows=2, free=T,connect="unique.pairs", values=acrossdS,labels=dz, lbound=-1,ubound=1.5), mxData( observed=corDZ, type="cor" ,numObs=365)), mxAlgebra( MZ.objective + DZ.objective, name="min2sumll" ), mxAlgebraObjective("min2sumll") ) ModelSAT$DZ@matrices$S@values <- ModelSAT$MZ@matrices$S@values + diag(c(0,rep(.7,11)),nrow=12) # since I was getting non-positive definite starting value matrices ModelSAT$MZ@matrices$S@values <- ModelSAT$MZ@matrices$S@values + diag(c(0,rep(.7,11)),nrow=12) FitSAT <- mxRun(ModelSAT) ### I also tried splitting up the covariance and variance elements with just as little luck within <- paste("a",1:10,sep="") withinF <- paste("b",1:6,sep="") mz <- paste("mz",1:15,sep="") dz <- paste("dz",1:15,sep="") acrossmS <- upperTriangle(t(corMZ[varsA,varsB]),diag=T) acrossdS <- upperTriangle(t(corDZ[varsA,varsB]),diag=T) ModelSAT1 <- mxModel("SAT", mxModel("MZ", type="RAM", manifestVars=vars, mxPath(from=varsA, to=varsA, arrows=2, free=T,connect="unique.bivariate", values=.5,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsB, to=varsB, arrows=2, free=T,connect="unique.bivariate", values=.5,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsA, arrows=2, free=T,connect="unique.pairs", values=.25,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsB, arrows=2, free=T,connect="unique.pairs", values=.25,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsFAM, arrows=2, free=T,connect="unique.bivariate", values=.5,labels=c("c1"), lbound=-1,ubound=1.5), mxPath(from=varsA, to=varsA, arrows=2, free=F,connect="single", values=1,labels=c("e1","e2","e3","e4","e5"), lbound=-1,ubound=1.5), mxPath(from=varsB, to=varsB, arrows=2, free=F,connect="single", values=1,labels=c("e1","e2","e3","e4","e5"), lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsFAM, arrows=2, free=F,connect="single", values=1,labels=c("c2","c3"), lbound=-1,ubound=1.5), mxPath(from=varsA, to=varsB, arrows=2, free=T,connect="unique.pairs", values=acrossmS,labels=mz, lbound=-1,ubound=1.5), mxData( observed=corMZ, type="cov" ,numObs=203)), mxModel("DZ", type="RAM", manifestVars=vars, mxPath(from=varsA, to=varsA, arrows=2, free=T,connect="unique.bivariate", values=.5,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsB, to=varsB, arrows=2, free=T,connect="unique.bivariate", values=.5,labels=within, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsA, arrows=2, free=T,connect="unique.pairs", values=.25,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsB, arrows=2, free=T,connect="unique.pairs", values=.25,labels=withinF, lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsFAM, arrows=2, free=T,connect="unique.bivariate", values=.5,labels=c("c1"), lbound=-1,ubound=1.5), mxPath(from=varsA, to=varsA, arrows=2, free=F,connect="single", values=1,labels=c("e1","e2","e3","e4","e5"), lbound=-1,ubound=1.5), mxPath(from=varsB, to=varsB, arrows=2, free=F,connect="single", values=1,labels=c("e1","e2","e3","e4","e5"), lbound=-1,ubound=1.5), mxPath(from=varsFAM, to=varsFAM, arrows=2, free=F,connect="single", values=1,labels=c("c2","c3"), lbound=-1,ubound=1.5), mxPath(from=varsA, to=varsB, arrows=2, free=T,connect="unique.pairs", values=acrossmS,labels=mz, lbound=-1,ubound=1.5), mxData( observed=corDZ, type="cov" ,numObs=365)), mxAlgebra( MZ.objective + DZ.objective, name="min2sumll" ), mxAlgebraObjective("min2sumll") ) FitSAT1 <- mxRun(ModelSAT1)