Login

News

Expected covariance matrix is not positive-definite in data row... at iteration...

Posted on
No user picture. dadrivr Joined: 01/19/2010
Forums
OpenMx Error Messages
Hi guys,

I'm trying to fit a saturated model where the variable, 'manifests', includes all of the variables in the model. The non-saturated structural model runs fine, but I get an error when I fit the saturated model:

Error: The job for model 'Saturated Model2' exited abnormally with the error message: Expected covariance matrix is not positive-definite in data row 79 at iteration 0.

I have tried fiddling with the starting values, but to no avail. By the way, I know that all=TRUE is being deprecated, and I have written the same code in a longer form without all=TRUE, but this was shorter, so I thought it'd be easier to post. Below is the model:

saturated_model <- mxModel(paste("Saturated Model"),
type="RAM",
manifestVars=manifests,

#Variances and Covariances
mxPath(
from=manifests,
arrows=2,
free=TRUE,
all=TRUE,
values=.5,

#Means
mxPath(
from="one",
to=manifests,
arrows=1,
free=TRUE,

mxData(observed=modeldata, type="raw")
)
saturatedModel_run <- mxRun(saturated_model)

Any help would be greatly appreciated. Thanks!

Replied on Tue, 09/20/2011 - 11:22
Picture of user. Ryne Joined: 07/31/2009

This is a starting values problem. No matter what constant value you pick for the single "variances and covariance" path, your expected covariance matrix will not be positive definite because all variables will be perfectly correlated. Instead, make one set of paths (all=TRUE) for the covariances with one constant, then make a second set of paths for just variances (all=FALSE) and give that set of paths a different, higher constant. Whenever you get this error at iteration 0 (the OpenMx backend is in C, which starts counting from zero), the problem exists in your initial starting values; higher iterations means that the optimization led to a bad expected covariance.

#Covariances
mxPath(
from=manifests,
arrows=2,
free=TRUE,
all=TRUE,
values=.5),

#Variances
mxPath(
from=manifests,
arrows=2,
free=TRUE,
all=FALSE,
values=1),

Log in or register to post comments
Replied on Tue, 09/20/2011 - 14:33
No user picture. dadrivr Joined: 01/19/2010

Great, that solved the problem. Do you know if there will be an alternative to "all=TRUE" for quick model specification in the next version of OpenMx? How should I code these models in the interim (i.e., omitting 'all=TRUE')? Do I need to code each set of paths separately? Thanks!
Log in or register to post comments
Replied on Sun, 01/29/2012 - 05:49
Picture of user. tbates Joined: 07/31/2009

"all=" has been replaced with "connect="

This is MUCH easier, more powerful, and easier to understand.

You can now connect from and to in the following ways:

"all.pairs": All possible combinations amongst from and to
"all.bivariate": As above, but excluding combinations from an parameter back to itself (covariances)

"unique.pairs": All possible combinations excluding mirror images (reversed from and to)
"unique.bivariate". "from" != "to" && no mirror images

"single": Just pairs up the from and to elements, as item 1 of "from" to item 1 of "to", etc.

Log in or register to post comments
Replied on Thu, 04/05/2012 - 02:55
No user picture. prast Joined: 08/18/2010

Hi, I have code which used to run with the all=TRUE argument but now, somehow, the exact same model does not run anymore.
I changed the all=TRUE to connect='all.pairs' and arrows=2 to arrows=1 (see old values behind the comment #). I would think that this should be the equivalent:

mxPath(
from=c("L1","S1","L2","S2"), # bivariate GCM with correlated intercepts and slopes
arrows=1,
#arrows=2,
connect='all.pairs',
#all=TRUE,
free=TRUE,
values=c(SigmaDep),
labels=c("varL1", "cL1S1", "cL1L2", "cL1S2",
"cL1S1", "varS1", "cS1L2", "cS1S2",
"cL1L2", "cS1L2", "varL2", "cL2S2",
"cL1S2", "cS1S2", "cL2S2", "varS2"),
lbound=c( 0, NA, NA, NA,
NA, 0, NA, NA,
NA, NA, 0, NA,
NA, NA, NA, 0)
),

After running the model I get 'The job for model 'bivgcm' exited abnormally with the error message: Expected covariance matrix is not positive-definite in data row 122 at major iteration 0 (minor iteration 1).'
Where am I going wrong?
Somehow the labels and values now showp up in the A matrix instead of in the S matrix. Is this intended OpenMx 1.2.?

Log in or register to post comments
Replied on Thu, 04/05/2012 - 10:22
Picture of user. tbates Joined: 07/31/2009

The problem (and the change from S (ymmetric) to A (symmetric) is because you have changed your paths from 2 to 1 arrow. You don't need (or want) to do that.

To replicate what all did by default but without the illegal mirror-pairs, just say:

connect = "unique.pairs"

To link each source to everything except itself, and exclude equivalent paths, say:
connect = "unique.bivariate"

The full option list is:

connect = c("single", "all.pairs", "unique.pairs", "all.bivariate", "unique.bivariate"),

Log in or register to post comments
Replied on Thu, 04/05/2012 - 09:22
Picture of user. mspiegel Joined: 07/31/2009

We are deprecating the ability of mxPath() to generate redundant specifications of the same path. This is why mxPath(c('a','b'), connect = 'all.pairs', arrows = 2) is not permitted, because the path from 'b' to 'a' is a redundant specification of the path from 'a' to 'b'. You can rewrite your path specification using connect = "unique.pairs" and arrows = 2:

mxPath(from = c("L1","S1","L2","S2"), 
    arrows = 2,
    connect = "unique.pairs",
    free = TRUE,
    labels = c("varL1", "cL1S1", "cL1L2", "cL1S2",
        "varS1", "cS1L2", "cS1S2", "varL2", "cL2S2",
        "varS2"))

You'll need to determine how to re-write your 'values' and 'lbound' arguments. Another alternative is to split this into two constructions, using "unique.bivariate" and "single". Sometimes it is easier to specify the 'values' and other arguments with this alternative construction:

mxPath(from = c("L1","S1","L2","S2"),
arrows = 2, connect = "unique.bivariate",
free = TRUE, labels = c("cL1S1", "cL1L2", "cL1S2",
"cS1L2", "cS1S2", "cL2S2"))

mxPath(from = c("L1","S1","L2","S2"),
arrows = 2, connect = "single",
free = TRUE, labels = c("varL1", "varS1", "varL2", "varS2"))
Log in or register to post comments
Replied on Sat, 01/28/2012 - 19:00
No user picture. brauer Joined: 01/28/2012

Hi,

I have to admit that I don't understand your answer. I used the script below to test a fully saturated model with four observed variables. I gave the paths a starting value of .1, the covariance a starting value of .5, and the variances a starting value of 1. The model runs but it produces an output that is different form the one I get with EQS or AMOS (different parameter estimates). Where is the problem? I tried to take out the mean structure or to rescale the variables so that they have more similar variances but this didn't help either. Any ideas?

-- Markus

data_apgar <- read.spss("data_apgar.sav", to.data.frame = T)
manifests <- c ("apgar", "gestat", "wgtgain", "smokes")
apgar_model1 <- mxModel ("Apgar Model 1",
type="RAM",
manifestVars = manifests,
mxPath(from="smokes", to="gestat", values=.1, label="b"),
mxPath(from="smokes", to="apgar", values=.1, label="c"),
mxPath(from="wgtgain", to="gestat", values=.1, label="d"),
mxPath(from="wgtgain", to="apgar", values=.1, label="e"),
mxPath(from="gestat", to="apgar", values=.1, label="f"),
mxPath(from=manifests, arrows=2, free=TRUE, values=1),
mxPath(from="smokes", to="wgtgain", arrows=2, values=.5, label="a"),
mxPath(from="one", to=manifests, arrows=1, free=TRUE, values=1),
mxData(data_apgar, type="raw")
)
apgar1out <- mxRun(apgar_model1)
apgar1out@output
summary(apgar1out)

Log in or register to post comments
Replied on Mon, 01/30/2012 - 15:37
No user picture. brauer Joined: 01/28/2012

Hi again,

Your help is greatly appreciated. I am not sure I know how to read the output. I pasted the output in a word document (see attached doc). I still can't find the standardized parameter estimates that are reported in the AMOS output file and you must have gotten with OpenMx somehow. For example, the path from "smokes" to "wgtgain" is -.461. Where did you find this estimate in the output? I searched the document for "-.46" (and -0.46) but didn't find it.

Then there are a number of things that are not in the output but that I probably have to ask for explicitly with the OpenMx script: The significance level of the paths, the standardized paths, the correlations between exogenous variables, the sample covariance matrix, the model-implied covariance matrix, the residual covariance matrix, the squared multiple correlations (how much variance my model explains in each of the endogenous variables). If you can tell me what I have to do to get the relevant output, that would be great. If not can you direct me to the appropriate documentation? Thanks a million in advance.

Log in or register to post comments
Replied on Mon, 01/30/2012 - 17:15
Picture of user. tbates Joined: 07/31/2009

Yes, you're right: the OpenMx summary gives a range of data, but you can, and should "ask for more".

Significance levels of paths are best tested with CIs (which will be included in the summary when requested in the model).
Standardization, I attach a script for (written by Ryne originally?)

Many things, though, are standard R: so sample size and count of manifests are just
dim(yourData)

Correlations and covariances, again you have R
cor(yourData)
cov(yourData)

etc.

Log in or register to post comments
Replied on Tue, 01/31/2012 - 10:02
No user picture. brauer Joined: 01/28/2012

Hi, thank you for your help. I looked at the description of MxPath and MxModel but couldn't figure out how to ask for confidence intervals in the model. Can you tell me how?

Standardization: Is there really no easier way than to write a 75-line script? My students are going to flip out if I show them the script you sent me! Getting standardized estimates seems such a frequent thing to do that I can't believe that there is an easier way to do it.

I still didn't entirely understand how to get the model-implied covariance matrix (I assume it's EC1 <- mxEval(S, apgar1out)) and the residual covariance matrix.

Log in or register to post comments
Replied on Tue, 01/31/2012 - 12:22
Picture of user. tbates Joined: 07/31/2009

Hi brauer,

To compute confidence intervals, you add mxCI() statements to the model, then call mxRun with intervals=TRUE.

mxCI statements are simply a list of labels, or matrices, or matrix cell addresses to caluculate Cis on.

So, from memory, you could say "mxCI( c('A') )" to get CIs on all the asymmetric paths, or mxCI(c('a','b','c')) to get CIs on your labelled paths.

The standardization function is lengthy. Because RAM is one arbitrary type of model that OpenMx can solve, it does not come packaged with all the functions that a person focussed on RAM-specified path models would use and want.

The idea is that you can build a small library of helper functions to automate the tasks you need repeatedly, hopefully sharing those here, like the standardization function.

Each group would then package up the (relatively few) functions they need, and source() those at the top of each script.

So if you drop the stand function and any others you make up for teaching into a file, that is instantly available to all your students as if it was built-in, but customisable to your exact needs and style.

I agree that a help file from the "mothership" as it were is a good idea, and I think that will emerge as the user base grows and settles on some standard shared requirements. Your questions help with that, so thank you!

Hope this helps,
tim

Log in or register to post comments
Replied on Thu, 02/02/2012 - 14:57
No user picture. brauer Joined: 01/28/2012

Hi again,

How did you generate the graphic representation of the model in the file you sent me ("samizdat.pdf")? Is there a way to get R to draw the graph for me? With the path coefficients? Thank you!

Log in or register to post comments
Replied on Tue, 10/08/2013 - 17:26
Picture of user. yenni Joined: 10/08/2013

Hi guys,

I'm trying to fit a saturated model. The non-saturated structural model runs fine, but I get an error when I fit the saturated model: Error: The job for model 'Saturated Model' exited abnormally with the error message: Expected covariance matrix is not positive-definite in data row 1041 at iteration 0(minor iteration 1).

I have tried fiddling with the starting values and did several suggestions from this forum, but it still failed,so I thought it'd be easier to post.
I post you the data set and my saturated file. I hope someone can help me. Thanks.

Log in or register to post comments
Replied on Tue, 10/08/2013 - 17:35
Picture of user. Ryne Joined: 07/31/2009

This isn't a saturated model. This is a common factor model with no residual variance terms. As all 50-something manifest variables are linearly dependent on the 9 or so latent variables, your model is not positive definite. Add residual variance terms for the manifest variables (the diagonal of the S matrix) and the model will be identified.

Why are you using this as a saturated model?

Log in or register to post comments
Replied on Tue, 10/08/2013 - 18:45
Picture of user. yenni Joined: 10/08/2013

Thanks for your response. I'm trying on it now.
Actually, I'm trying to compare OpenMx results and LISREL results according to the paper of Statistica Neerlandica(Toharudin, Oud, Billiet 2007).
The structural model runs fine, although some parameters was not same with LISREL results.
My problem is I can not calculate chi square, so I try it to calculate Saturated Likelihood.

I have also tried with syntax below to calculate it, but it keeps running forever.
###compute Likelihood of saturated model and compare with likelihood of postulated model
Sat.Model <- function(data1){ # data1 are the data WITHOUT the time intervals at the end
SatModel <- mxModel("SatModel",
mxMatrix(type="Full", nrow=1, ncol=ncol(data1), free=TRUE, values=sapply(data1, mean, na.rm=T), name="expMean"),
mxMatrix(type="Symm", nrow=ncol(data1), ncol=ncol(data1), free=TRUE, values=var(data1, na.rm=T), name="expCov"),
mxData(observed=data1, type="raw"),
mxFIMLObjective(covariance="expCov", means="expMean", dimnames=names(data1))
)
SatModelFit<-mxRun(SatModel)
summary_SatFit<-summary(SatModelFit)
return(c(summary_SatFit$Minus2LogLikelihood, summary_SatFit$degreesOfFreedom))
}

chi_square <- function(data1, summary_fit)
{
temp <- Sat.Model(data1)
LLsat <- temp[1]
df_sat <- temp[2]

summary_fit$SaturatedLikelihood <- temp[1]
summary_fit$saturatedDoF <- temp[2]
summary_fit$Chi <- summary_fit$Minus2LogLikelihood - LLsat
chi_df <- summary_fit$degreesOfFreedom - df_sat
summary_fit$p <- (1-pchisq(summary_fit$Chi, chi_df))
summary_fit$RMSEA <- sqrt(((summary_fit$Chi/chi_df)- 1)/summary_fit$numObs )
return(summary_fit)
}

at the same time I'm also running your suggestion, but it keeps still running. I'm afraid that the result will be same with the previous one.

Log in or register to post comments
Replied on Tue, 10/08/2013 - 19:18
Picture of user. mhunter Joined: 07/31/2009

Try something like the following:


###compute Likelihood of saturated model and compare with likelihood of postulated model
Sat.Model <- function(data1){ # data1 are the data WITHOUT the time intervals at the end
SatModel <- mxModel("SatModel",
mxMatrix(type="Full", nrow=1, ncol=ncol(data1), free=TRUE, values=sapply(data1, mean, na.rm=T), name="expMean"),
mxMatrix(type="Symm", nrow=ncol(data1), ncol=ncol(data1), free=TRUE, values=var(data1, na.rm=T), name="expCov"),
mxData(observed=data1, type="raw"),
mxFIMLObjective(covariance="expCov", means="expMean", dimnames=names(data1))
)
return(SatModel)
}
mySaturatedModel <- Sat.Model(myDataToFit)
mySaturatedFit <- mxRun(mySaturatedModel)

myComparisonModel <- mxModel(
#...
)

myComparisonFit <- mxRun(myComparisonModel)

summary(myComparisonFit, SaturatedLikelihood=mySaturatedFit)

If that does not work, then let us know where in the process it failed and we should be able to fix it!

Log in or register to post comments
Replied on Tue, 10/08/2013 - 17:27
Picture of user. Ryne Joined: 07/31/2009

So this is a starting values problem, and one caused by your use of the "all" argument. You're creating all of your variances and covariances with the same function call, so whatever starting value you assign to the variances gets assigned to the covariances as well. That means that you've assigned all variables to be perfectly correlated.

I'll also point out that the "connect" argument is much more flexible than the "all" argument. connect="all.pairs" does what "all=TRUE" does, creating all possible paths for a set of variables. For two-headed arrows, that's actually a pain, as it creates a with b AND b with a. You have four options for connect, which will serve your needs well.

all.pairs does what all=TRUE does.
unique.pairs will do all variables with all variables EXCEPT repeat the redundant paths. a with b will be created, but b with a won't.
all.bivariate and unique.bivariate do the same thing, but leave out the variance terms (a with a is dropped).

So what do you do? Use two mxPath statements. In the first, leave out the all and connect arguments and create the variance terms. In the second, use connect=unique.bivariate to make all of the covariances. Make sure the start values for the covariances are lower than the variances and you're all set.

Log in or register to post comments
  Login

News