OpenMx from 0.1.3-776 to 0.1.4-827, what changed?

19 replies [Last post]
pdeboeck's picture
Offline
Joined: 08/04/2009

Hi everyone,

Here's an interesting problem I am running into. I was working on improving convergence rates on a script that was giving me some problems. I finally got better convergence on my laptop after switching to FIML and changing the scaling of coefficients. (I think the variances of x and derivatives I was estimating were differing by several orders of magnitude.) Then it was time to get things working on Xgrid. After a while I noticed that one of my Xgrid computers (my desktop) was almost never converging, while my other computers were converging almost all the time.

There were two differences that I could find: my desktop had R 2.9.1 (other computer 2.9.2) and my desktop had a new version of OpenMx (other computer was updated 2 weeks ago). I updated R on the desktop, no improvement. I updated OpenMx on the other machine --- almost none of the models converge now on both machines.

Something has changed with OpenMx in the new version. Is there something new that I might not be specifying in a function which was not previously required? Could there be something different with mxFIMLObjective? Other ideas?

pdeboeck's picture
Offline
Joined: 08/04/2009
Question unrelated to the

Question unrelated to the previous topic, but related to changes in the new version:

I am running a large number of models using a for loop. Before (v0.1.3) when running multiple models, R didn't tend to give errors that would cause it halt execution. Consequently, I could run through hundreds of models without problem. (Naturally I am storing the status codes so that I can check that models converged.) Now I get the following: "Error: The job for model 'LDE_Model_1' exited abnormally with the error message: Non-positive-definite at free parameters: 0.793 -1.439 0.072 0.054...." This error terminates the current loop, canceling the fitting of numerous additional models. Is there any way for OpenMx to register this in the status codes, rather than halting execution of the program?

mspiegel's picture
Offline
Joined: 07/31/2009
You could use result <-

You could use result <- try(mxRun(model)). Keep in mind that if an error occurs, you still can't get the result of the optimization under the current system. If an error occurs, then the value stored by result is an invisible object of class '"try-error"' containing the error message.

pdeboeck's picture
Offline
Joined: 08/04/2009
Somehow I never ran into

Somehow I never ran into try(). I'll give a try. Thanks!

pdeboeck's picture
Offline
Joined: 08/04/2009
Here's some sample code.

Here's some sample code. Sorry about the .R extensions on everything, the forum would not let me use .dat. Run Test.R, it will load information from the other two files.

server: 2 code 4's, 23 code 6's.
desktop 2: 2 code 4's, 23 code 6's.
laptop: 23 code 0's, 2 code 1's.
desktop 1: 23 code 0's, 2 code 1's.

All are running R 2.9.2. The server and desktop 2 are using v0.1.4 of OpenMx. The laptop and desktop 1 are using v0.1.3 of OpenMx.

AttachmentSize
Test.R 2.18 KB
loadings.R 203 bytes
Data.R 50.2 KB
tbrick's picture
Offline
Joined: 07/31/2009
Did we change the effects of

Did we change the effects of byrow=TRUE?

Try removing byrow=TRUE from the line about matrix L (line 42, I think).

This fixed it on my machine.

pdeboeck's picture
Offline
Joined: 08/04/2009
Looks good initially --- I'll

Looks good initially --- I'll run a larger test tomorrow.

For those who might be interested: The values of matrix L were set by a E x 3 matrix of the same name. It looks like the L matrix was correctly specified in v0.1.3 regardless of whether the "byrow" statement was included. (Even setting it to FALSE resulted in the same matrix being produced.) In v0.1.4 the inclusion/removal of the byrow statement makes a difference.

pdeboeck's picture
Offline
Joined: 08/04/2009
The test went well, aside

The test went well, aside from not getting to it until Friday. The change in the L matrix more than doubled the number of converging models. Thanks!

mspiegel's picture
Offline
Joined: 07/31/2009
Poop. There is a bug in

Poop. There is a bug in 0.1.4-827 and the bug is that when 'values', 'lbound', or 'ubound' are specified as matrix arguments to the mxMatrix() function AND when byrow = TRUE then the matrices are read by columns. Let me explain how on earth this happened. A few weeks ago I noticed one of our scripts had the call mxMatrix(values = 1:9) and this causing OpenMx to report an error because 1:9 are of integer type and values must be of float type. So I added the line to mxMatrix: values <- as.numeric(values). And this line is AFTER I test that is.numeric(values) returns TRUE. * Quirk #1 of the example: is.numeric() returns true on integers or floating point numbers, and as.numeric() always returns floating point numbers. What I did not realise is that as.numeric(values) will convert a matrix into a vector (!?!) Furthermore it converts the matrix into a vector in column-major order, such that when byrow = TRUE the column-major vector is read in row-major order and the badness happens.

The fix has been checked into subversion, along with a test case to make sure this doesn't happen again. We'll release 0.1.5 with this fix on Friday.

Steve's picture
Offline
Joined: 07/30/2009
Glad to hear this one was

Glad to hear this one was found. Nice catch.

pdeboeck's picture
Offline
Joined: 08/04/2009
Ah - that clarifies my

Ah - that clarifies my observation (written below).

I never did like that as.numeric() would change matrices into vectors.

tbates's picture
Offline
Joined: 07/31/2009
so that's bad for me here

so that's bad for me here too

SVN Mx, "R version 2.9.2 (2009-08-24)"

Running LDE_Model_1 
Rep:  1 	 Wed Sep 23 15:28:06 2009 	 Status:  6 
Running LDE_Model_1 
Rep:  2 	 Wed Sep 23 15:28:11 2009 	 Status:  6 
Warning messages:
1: In model 'LDE_Model_1' NPSOL returned a non-zero status code 6. model does not satisfy the first-order optimality conditions to the required accuracy, and no improved point for the merit function could be found during the final linesearch (Mx status RED) 
2: In model 'LDE_Model_1' NPSOL returned a non-zero status code 6. model does not satisfy the first-order optimality conditions to the required accuracy, and no improved point for the merit function could be found during the final linesearch (Mx status RED) 

pdeboeck's picture
Offline
Joined: 08/04/2009
Earlier you mentioned

Earlier you mentioned downgrading. Have you tried that? Does it fix the problem?

pdeboeck's picture
Offline
Joined: 08/04/2009
If the symptoms help: The

If the symptoms help: The laptop (OpenMx 0.1.3) with the weaker processor has no problem outpacing the desktop (0.1.4). I just ran 10 models, the desktop produced 10 code-reds, the laptop produces 10 zeros (i.e. all was fine).

tbates's picture
Offline
Joined: 07/31/2009
So you downgraded the laptop

So you downgraded the laptop again (you mentioned you tried it on the new version and got failures to converge)?

What happens when you downgrade the desktop?

Do you have the script to share with a url for the data?

pdeboeck's picture
Offline
Joined: 08/04/2009
This is based on the script

This is based on the script in the following forum:

http://openmx.psyc.virginia.edu/thread/144

I haven't tried that script, but it should produce similar problems.

tbates's picture
Offline
Joined: 07/31/2009
hi, Running fine here:

hi,
Running fine here: returns in about 1 second, no errors with openmx from SVN and R 2.9.1. Updated to 2.9.2, and still molto bene.

ldeModel1Fit <- mxRun(ldeModel1)
Running LDE_Model_1
name matrix row col parameter estimate error estimate
1 A 3 1 0.041294837 20.0094740
2 A 3 2 -1.206272186 4.8355851
3 S 1 1 0.134677338 3.2759996
4 S 1 2 -0.231808468 6.1104767
5 S 2 2 0.000001000 4.4628945
6 S 3 3 0.000001000 3.2374304
7 U 1 1 0.018200208 0.1093465
8 U 2 2 0.003239532 0.1314263
9 U 3 3 0.016667088 6.8067205
10 U 4 4 0.129030656 0.7500015
11 U 5 5 0.149651825 3.4870282
12 U 6 6 0.116827792 2.2173529
13 U 7 7 0.145345958 1.7590130

Observed statistics: 29
Estimated parameters: 13
Degrees of freedom: 16
-2 log likelihood: -1204.045
Saturated -2 log likelihood: -1765.218
Chi-Square: 561.1724
p: 3.871965e-109
AIC (Mx): 529.1724
BIC (Mx): 244.2398
adjusted BIC:
RMSEA: 0.6020645

running under 2.9.2
Observed statistics: 29
Estimated parameters: 13
Degrees of freedom: 16
-2 log likelihood: -1204.045
Saturated -2 log likelihood: -1765.218
Chi-Square: 561.1724
p: 3.871965e-109
AIC (Mx): 529.1724
BIC (Mx): 244.2398
RMSEA: 0.6020645

pdeboeck's picture
Offline
Joined: 08/04/2009
I'll post a script that

I'll post a script that recreates the problem later today. (First I'm off to prep & teach a class.)

pdeboeck's picture
Offline
Joined: 08/04/2009
I hadn't gotten around to

I hadn't gotten around to updating the laptop. I've been using source() to get the updates, so I'm not sure how I would go about downgrading the desktop.

The code includes a lot of other pieces, but here are the pieces related to the mxModel:

ldeModel <- mxModel("LDE_Model_1",
mxMatrix("Full", 3, 3,
values=c( 0, 0, 0,
0, 0, 0,
-.01,0, 0),
free=c(FALSE,FALSE,FALSE,
FALSE,FALSE,FALSE,
TRUE, TRUE,FALSE),
name="A", byrow=TRUE),
mxMatrix("Symm", 3, 3,
values=c(.1,0, 0,
0,.1, 0,
0, 0,.1),
free=c( TRUE, TRUE,FALSE,
TRUE, TRUE,FALSE,
FALSE,FALSE, TRUE),
lbound=c(0,NA, NA,
NA,0, NA,
NA,NA,0),
name="S", byrow=TRUE),
mxMatrix("Iden", 3, name="I"),
mxFIMLObjective("R","Means")
)

RunModel <- function(toAnalyze,loadings,manifestVars,type) {
if(type=="Orig") L <- loadings$LDE.Orig
if(type=="GOLD") L <- loadings$LDE.GOLD
len <- length(manifestVars)

fit <- mxRun(mxModel(ldeModel,toAnalyze, mxMatrix("Full",free=FALSE,values=L,name="L",byrow=TRUE), mxMatrix("Diag", len, len, values=.05,free=TRUE, name="U", lbound=0), mxMatrix("Full", free=TRUE,nrow=1, ncol=len,values=0, name="Means", dimnames=list(NULL,manifestVars)), mxAlgebra(L %*% solve(I-A) %*% S %*% t(solve(I-A)) %*% t(L) + U, name="R", dimnames = list(manifestVars,manifestVars))))

Out <- c( fit@matrices$A@values[3,1:2], #eta/zeta
diag(fit@matrices$S@values), #derivative var
fit@matrices$S@values[2,1], #derivative cov
#fit@output$SaturatedLikelihood, #saturated likelihood
#fit@output$Minus2LogLikelihood, #-2LL
fit@output$status[[1]]) #Error Codes

return(Out)
}

pdeboeck's picture
Offline
Joined: 08/04/2009
Oh, and

Oh, and this:

mxData(Embedded,type="raw", numObs=dim(Embedded)[1])