Hi all,
I’m trying to run a saturated model with restrictions in a multivarite analysis with 6 continuous variables.
In this saturated model with restrictions (see below) I had the following problem (that I did not have in the saturated model without restrictions):
Error: The job for model 'Con' exited abnormally with the error message: Expected covariance matrix is not positive-definite in data row 431 at major iteration 0 (minor iteration 1).
However, when I use the same script but instead of 6 with 5 variables or less, the script runs perfectly. I don’t know what it’s happing… I think it’s something related to the starting values, but I don’t know how I can solve this.
Any help would be very helpful.
Many thanks in advance.
#
Restrictions: means and variances equated across birth-order & zygosity groups;
One set of Cross-trait cor, symmetric Cross-trait Cross-twin correlations (MZ and DZ)
Create start values
To avoid starting the optimization at the solution add some random noise
jiggle <-rnorm((nv(nv+1)/2), mean = 0, sd=.1) #calculo jiggle necesario segun estos tamaños
Stmean <-colMeans(Data,na.rm=TRUE)
Stsd <-sd(Data,na.rm=TRUE)+jiggle[1:nv]
StWithinperson <-vechs(cor(Data,use="complete")) +jiggle[1:(nv(nv-1)/2)]
StBetweenMZ <-vech(cov(mzData[,1:nv],mzData[,(nv+1):(2nv)],use="complete")) +jiggle
StBetweenDZ <-vech(cov(dzData[,1:nv],dzData[,(nv+1):(2nv)],use="complete")) +jiggle
Create Labels for Column and Diagonal Matrices
mLabs <- paste("m",1:nv,sep="")
sdLabs <- paste("sd",1:nv,sep="")
Create Lables for a Correlation Matrix
rphLabs <- paste("r",1:ncor,sep="")
Create Labels for Lower Triangular Matrices
MZbLabs <- paste("mz", do.call(c, sapply(seq(1, nv), function(x){ paste(x:nv, x,sep="") })), sep="")
DZbLabs <- paste("dz", do.call(c, sapply(seq(1, nv), function(x){ paste(x:nv, x,sep="") })), sep="")
Specify Matrices to hold parameter estimates for MZ and DZ data
elements for one overall Means
meansG <-mxMatrix("Full", 1, ntv, free = TRUE, values = Stmean, labels=c(mLabs,mLabs), name="ExpMeans")
elements for the SD
sdZ <-mxMatrix("Zero", nv, nv, free=F, name="padding")
sdG <-mxMatrix("Diag", nv, nv, free = TRUE, values = Stsd, labels=c(sdLabs), name="SD")
sdT <-mxAlgebra(rbind(cbind(SD,padding), cbind(padding, SD)), name="SDtwin")
elements for the correlations
Rph <-mxMatrix("Stand", nv, nv, free = TRUE, values = StWithinperson, labels=rphLabs, name="within")
MZb <-mxMatrix("Symm", nv, nv, free = TRUE, values = StBetweenMZ, labels=MZbLabs, name="BetweenMZ")
DZb <-mxMatrix("Symm", nv, nv, free = TRUE, values = StBetweenDZ, labels=DZbLabs, name="BetweenDZ")
corMZ <-mxAlgebra(rbind(cbind(within,BetweenMZ), cbind(BetweenMZ, within)), name="RMZ")
corDZ <-mxAlgebra(rbind(cbind(within,BetweenDZ), cbind(BetweenDZ, within)), name="RDZ")
generate expected covariance matrices
covMZ <-mxAlgebra(SDtwin %% RMZ %% t(SDtwin), name="ExpCovMZ")
covDZ <-mxAlgebra(SDtwin %% RDZ %% t(SDtwin), name="ExpCovDZ")
Data objects for Multiple Groups
dataMZ <- mxData( observed=mzData, type="raw" )
dataDZ <- mxData( observed=dzData, type="raw" )
Objective objects for Multiple Groups
objMZ <- mxFIMLObjective( covariance="ExpCovMZ", means="ExpMeans", dimnames=selVars )
objDZ <- mxFIMLObjective( covariance="ExpCovDZ", means="ExpMeans", dimnames=selVars )
Combine Groups
pars <- list(meansG, sdZ, sdG, sdT, Rph)
modelMZ <- mxModel(pars, MZb, corMZ, covMZ, dataMZ, objMZ, name="MZ" )
modelDZ <- mxModel(pars, DZb, corDZ, covDZ, dataDZ, objDZ, name="DZ" )
minus2ll <- mxAlgebra(expression=MZ.objective + DZ.objective, name="m2LL" )
obj <- mxAlgebraObjective( "m2LL" )
ConCorModel <- mxModel("Con", pars, modelMZ, modelDZ, minus2ll, obj )
ConCorFit <- mxRun(ConCorModel)
ConCorSumm <- summary(ConCorFit)